


‘ 


‘ Ry. 


‘W S aN 
iby F ‘ hh Me , 









Institutional Archive of the Naval Postgraduate School 


Calhoun: The NPS Institutional Archive 
DSpace Repository 


Theses and Dissertations 1. Thesis and Dissertation Collection, all items 


1988 


Adaptive algorithms for two dimensional filtering. 


Wilstrup, Steven L. 


Monterey, California. Naval Postgraduate School 


http://ndl.handle.net/10945/22855 


Downloaded from NPS Archive: Calhoun 


| Calhoun is the Naval Postgraduate School's public access digital repository for 
e (8 D U DLEY research materials and institutional publications created by the NPS community. 
«ist | et Calhoun is named for Professor of Mathematics Guy K. Calhoun, NPS's first 
TT \ KNOX appointed -- and published — scholarly author. 
http://www.nps.edu/library 






LIBRARY 


Dudley Knox Library / Naval Postgraduate School 
411 Dyer Road / 1 University Circle 
Monterey, California USA 93943 


Oo Ao 
1 
o ' 
7 5 
‘ 
i} 
A 
iy * 
O 
O 
fo 
D 
nT 
O 
OD 
Ln 
' O 
' A 
. Oo 
A 
a A so «4 
A ' 
D 
fo 
pO 
i] 
Pu 
D 
Ly 
oa . 
iy 
« 
* 
' 
J 
fo rae 
rd 
’ 
sl 
’ 
' 
Oo 
O 
D 
A 
iy 
D 
' 
cy * 
oy 
f 
A 
Dy 
1 
. 
1 D 
oO 
0 
iO 
. . 
O 
fo 
' « 
O 
J 
* 
. o 














~ im > Li MAb be 4 a8 = — hall seey a, bat coadie SO Oth ot te been eed Path ehee ht ehh ee PT ee 
re es ae eae poe ce ar) ye Oy eo eed er a eae me fon rr ppp teri e terion, beta tee et hee 
= A Tie eee t So he a ae Were ah) P aati coal [iatiiaheethhd de ute ee et 
are : : Bo dl Sk Rey eo Paper ea rei eee 4 np ‘ 
A 7 on wa) er er a 1 ena OR Re Lr ead TO es oe re Ee ao nee a. fn Piel LIE Tt Tile eee PEP ey he ne 
ar a) n cot re ee Si as eee te re ee Soe SP eee Tae a hie er ee 
or) A Cr an) n re ae ag he oe Se oe?) ee ee > Poll Pree oe ieee ee ot et ey i we eee ee te tt 
5 4a aes # ee rs ee ee Pe es ie tS Pe ese ee einen PLT ee ey ee 
’ ’ a Cn eh ee eS he ee et en ee ee ee Sn be ee eT re ee er ei Pier + 
A Q@eaur i i i. i, a ek Oe ee ae Le 7 ee ae re te po Tey 
‘ . a) ‘ n a rn he Coe |e) he hI) MT tat ee nd et Saab ieie OL e on bee ed ee oe Te Pee eee) 


ee « ar Pr CM Oe te Poe ce te) Pe 


a ee i ery) ee ree SU ee ot) ese the hore ay eee ee eee ea Ede et tT ety ery 


©: ener Sf Se A Ra he re earn ae Ran 

je 0 Me at ee LAL Seen ert oe eee ie Leeper neice ete al tla dle hein rt 
Si a eT A Le ee ee Le ee ee 

ee ee eC Mr et ee CL kt Pe ne er ee 

e Ce Pt Tt oe © ae, ct ee) oo verry) Peo ee vie peasy 

U Phe a pCa ry | ae Ce ee ee ee ee Ve ee ee ed oe ee 

td = s88 q See OF Pee Te) Se © aed he Er i er eee Pe ee 








aw.cs « 
a Ti ee iy eo re ts er eee er Pry) er ary dh nineteen he Lie ee 


lied tie letedt eh Y Ms 
a Ce eRe MRL Loe Sa ee Or Oe Pee ee tee rb tehhel tart) hapten yn helen 





ee eh I Pe 
ST ihre. f La) 
PLease 
es a.m 


me | i Ler? de i oak Ty , 
ny a en ee ee sy 
D [a oe a a) ee 
“ot! omg @, 1h" o fee? 0% 4 LAWene 
. fr ee 2 ee hated ees heer pred ae 
an] rd o TL? Tt St Lae Ie id Fy Shad nh] Sa aire o 
ie ; toh aie o nines om Veh REG FON RAE QIN obt@cen eds apunetgnn i ekeditt Teta) cnt ee 
, a ey ee Ce oh SC CS Pon ECT ob? et. Ym Pl 
eT ¢* @ 4 Cet ee Pe ee tt bet TT) 

n D eo @ a8’ "0989 SOY Cee ay otey B40 ews orp 
Ke Ps Pardo pa + em ery ee 
¢ 6A @as of ee. 8 PY Th at ee Led ae # org Po Sed er es ett wh ee 

O 18 6 thes 0 8 BM Otegs Tate e fo Abdi OR | bed SE TRS OPTOMA 5 1 oa nde 
Ceo ‘eo © ee ore ee ee hk ee er Pi tel ee cos lt veer eC Ty nape as ES 







halt tol tT Te aa ee oe 
yore) Pe rete] any Tene Pea hepa tar a- bean 
eACOMOT BBR AS ein Re hide de ees ee Tee 

tee eT ae LO eel Peel TT pe re ay eee eee 

























an AX +. ’ Py sf ater ree We ee ee woe ees ee aS ee ee terre vy oe ee Perr Pret i 


' at . P'e@ee , 4e | ee 2 eT ee ee te Ly ie Yo) me eid 


Ca | . . i - t@eoee hh ee er Pe he Tred 


a o4* 6 "4 


tt be ng 8 ate a ty ote ie 8B, 8 Oat OFO rat rnit ams 69 OPP ts Qa e. og8 ogres dwaqe 4.5) 0-8 @ Rimek 


be dad bP Ld A 
See ee to ee rT YO CS Ee Mees Me ot Ie eT SE a eT ee pee caper Plana o> 
® Li ee ee ee ee a Cr Chor eee Mee iy ye ea any aioe ee on 
Oar. a, Ce ee roe Mey Le a ee eee ee pe erere AT apd ppl ah 
a i a LP ee ee er et) er ee ey bh Te ee 
aterece™ bahnadiil teintl antes bee Te eer? 
Cae oer ee bees eo re Pee ee er pre beaengl 


Pate hb ie | 
Pe as oT Col ee er tk Sree pear pry ean Hh haem 


Ch ee Oe i eS i ee Sh ee Ted fe ee Pe ete ee = ‘ae 
a 2 tof eghs oe. Ce ae ed wdouba-dgeane Pee a es Loddedtitel eam! PE etter ptibegeen 
th Ne Miia! IG aS hee a shel ee ED el ol led ats IRR ee th Lae hae eee Sin 

Le oe eee he ed ee te ey ees be diy fe la al be gate eg 
ee eed Tl) aie Che Ue eee ll Pee i eT ee eet Te 

Oy ak ee ey heed ete Sie Th ee te Pe ee | hoe ae 

eo ed iy ee ee tg Ce eee ey eee a 

Cec met t] o % r fo Were Cob st, Verte © oO, enfer 4h. Morne’ ds o%@ dete te ech LLL ear eet ee ren 

oe F » ue Oe ee ee | eke Prd # oP O00 BO prsH sm, tg SPX asad g bd dental eel ial enh eh Lal Dat Ln Eee en fe 

F PPL i ie kk ee es Pe oe i rey ee Rae rit ee ee Re CL fe en en ane 

* es ie te Peer ie te er Pe ie yi yee Sime YP ee eee hte Melati eT ET Te 

“ann ee smse unre ge a hy ih de bh OL eel et RT ott eee Pat eet dt b letetla hated te atlina hare Tt Te Pe eee 
Pr i ee eh Ue hd er eT retn Sere There yy dite tae tent re popes 
er Pera tat eye he moh Ot ee) Dt Pre ee Cr eee ba Se ee eae 





on re ere ee *a Alera on 


















of a nn | « 


















br el meek 
baa Mel hh Tpit 
b latina held dae he 
belt et) A ee ee Th, 


A Cede PorTeret yy et, ee eee her eh yk ere ry ete Pee eng Ree 
® rey y ste an meres Be Be = BOlt Hee Wen? + pert 3.%.1 & SOF oe, By Ce ee 
. ~ 






: as OL} ' a Ch 


on . > eet hoe Sb] 

Po ec To ee Le oY bth ie at eet he CP re ee e 
CP ee er Oe. ere Tt or mr alc Ter | kel tl eM lace ee ye ey aay ey" 
Crna Ty Ore Sree oe eS er). ee ih ok hee el el ee ee] aT. 
CT ac ee See ae « a re ee ry s erta fe eb Ad eC bee] bd te eee ee eT 
i? he Ce oC oe ro. creer Te Serre yar er et eae hibits etek Th Tete ee 

ee aA a eT CM Ye) Pc Ce ed Ce a ee betta eT ELE cea Lee nyt eee ae 

Py eee f i eh ee Li oe et ee ee Ole Tip hy Cet RRY Ct Pee) 

Pe Te ee ee TP ry fe wn Aw Otte © Me ty te epee de iied inte Tey ee ens 
eRe MS te FE Fee Bt Ot Det Sater t- jsaceea! eprsesotee et ee eo eee rt 
: Poot RS 1g8e w& Ce rir) wetdrre Swe eee, PL eee ee ee tt 
AMC aracces Coat FAC Ce ee ee ee Pt ee ee et tT Peery Che aes ee ee ee 
- a Cr ee CO i rT) ht Ie nee Th it dad, Coe | phe betel Th ee ee 
+ pends re] - t« ers el? eee eee en Tire er | Pee | a Re peepee eden taper seer eres 4 
CU sosergt 2 © Sasne be © Me" Sewers ce bh lettin late! beled tel attendee ae ee 
‘4a a O a ht i ached oad lhl debe id ot tee ee 
2 ee e's re i delet. Tek tiie et ee th 
SoA ae. i he Se Ae Pee et Tir ey Ce ee ee 
a i | Leah be eel be ete a tL UTE LRT Le ee ne eee 
ee andthe hn ell Steet net rent Co ten 
hhh ht el cee eas 
Ore oe Oe REE Re ome Nis Oe oe be spe patac aecRa Y be 
Rit ae bo T) es ale UT Tk ee ee Oe 
: gfe oe eR ey eee ee Sh ae ee ae the p25 bere panne acme hates oa tech aol. 
Pra Oe | ee tee et Pk ee ne rr nee rad . 
Ps yt erehogt 2 o & se tabe Op fe bee obe cece th tee ine Manila th ele LEME eink oe eer 
Pre i ee el ee i er od Seccaugrres wrest s heen bate t-te be bh he ee 
a oe) i he tT ot Vebered gen) bere cweecas: OteeX teres eet tt er 
Ptwrn Pc per eee Seo = oh aed eolnn-s Me nt 41S wroqsteatous ea geome 
He aa oe er crt <] oC bt Oe OT LS hy en Teas 
Pe Ve ie ee ee Pe er TT | ee en © Fen Het eretaans 
‘ ae ce er he Ch et Tht ye or De 
@ Sqtsm 2 oe oe - % re fe. o 
rs 2 eten tf ote 






















































ontiae 
ore ty ae eee 





















ee te. et te oe re 
Sl te he at et ee 

am Yt eerl Serr rye ee pee 
eben te ated te It rn Tee ra 
a 








5 a 
UR eee ee Te ae ape 
“18 toe Dear eer een Per ore 2 
_ Yi cet oe Peer er Ree Part oe eer ern kT | Ce en 
A fee ca ee ee ee Ce tS pdabhecke <1 | Cs Sere ee ae ee ne 
oa lied 2 Ten BEC i i hee ht LS el Dh tei het te et Ce tee 
Fi = 8 rn | Sak ales ry a oT es ok ihr ae he 1) - ot tom Eerie eS 
= rn ® or ee ee ae Te Cr a. ? oe Se ee ern #4) . inp we take Pr pele ie 
« « Pes ees e info ce ee 4 a eS Ce allied a seer vot @tee Qe“ ate oh Me ngh en ee nie ar 
Sire Car Pe CL Ye Se Pit ok a ea al ae ary ya 
ry ' oP a ei) ey) a i ar) Oe at) °7 tt) Pe i ere ie oe ee in ee 
Ce a) ory ‘ev dint 0-8 e bo a! WG ay 5 SS Sao we Ber tsare- 2-8 rie * 
o 5 ans eh Pt Le Le Lee kere es POR R oe amt Fe own ela 9 ie 
Pi e rer ic re a S rn So * he Peery TE ee Ce ae 
rs Ae al: ee ee | ee) oe te ey 
one Pa? es ; 
Pe ae a * a San teh ERX DR whet Meal ar) 
+ ope 6 & te - rT errs Ee ee es Por 2? fee 2 
eT SS Bela ot PO ee ee oe Se ee 
Ty ee ee Crit a ere ee ee ee ee ree Saeein 
2 4 ey os FSS Aa pat DE deal Ba 
\ rel \ ides Pe Berta Pa Perera a 
ae a « © Se & 28 = rie ci nes . teh. 7 2 si 
eee a eet ier yma oer tia tale bt Pry Pee ares trent 
hh) re a he wee Gat ab? 0b 9 BENGE WIE es 
Ww a Featast he ae 8s 08S © hE mw cece i ded Pur, atte TRS te 
hr ee Pee) Coe ee Ds OV Le Se ee en TO oe Ce ee 
ee Sa cia ee Merk LL) SL ek tek of ae © a + pe Se 
' oo ow FY ' J ee os f Py a be ee dy ot g¢* “.° 4 
rer sad Li a Pr St hae De, SS oe 
tees Pe 2 er PLLC ce) Co a <9 
Py . . o e . < @sa%, 
le a i Par 


Ct) 


a7 
ks br ba he ihteele 2 a 














































bd ee 


a i ae 





























































ies et Oe ec ee 
: : ar See eee 
ae th en A 
re : ote Sita a a a> 
oh a > : ete ee ie Pere te 
: 7 i et Pr ae aed 
4 e i « Be Poe Pe 
; : es ty 
ae 
+ eer 7 
°, a | 
$ 7 Lagubd bet Sell 
rw t ‘ A = a 
' - 7 a. ee 
; pA be pelts «pa 
. 7 7 a a ht he Ln ap 
‘ Pe ed ad cy 
= a) Hs tet Ad tr eee CS ee 
5 a A « a et i i cee t,t eee 
; ss - a tote eh te 
7 Pr Mien ee ods fe 
; hee whe To ee 
an 5 es ee Sere aa ates 
: - * ° Ss dn et tract. dang Ry 
; bad nt er ey ee ee 
* 4 | «RO en “5 A2esten, et . o ree Lae 
oe Tee inet date ste Ty a Ree 
a "Fie . Ma She Mate al al b 
aaa a Es OP Sal Pe Sore abo 
5 Pa Co er a he eS Pee oe 
s ' ih. ie a ied LF ah, i hier Shane Se eee 
es 22) s aigt ta oF ial a: Se ee 
. wt Css me 4 BO Se ee ee 
' ' ¢ 4 bh bel te Clee bl oe . Bea = 
a) ° n A Lit, flat lie beth eee ¥-3 i oe eae 
8 . Pa i Weg tt Bir OO seed ieee eet tate’ nt Scat ae 
Hy bhtetadetel thera, Me tela 
hg ee behind i ie tee Wee 
J 
aa st belied eh th hd Steal ded bel fied Ste TL 
. . Ih filles fhe cells Pood Jon) Dah od tee stele 
: eS wie eat! ok 9 cole tat ed Mt Pe 
a a A ent nD Soatat Ske Pee) Ph fas 
ary ‘ A . scaled pl P pee ake ert Ome a, e aueree ye 
er] be bella hal tdel Lok te Ll CT eld = ° 
oar] 5 28 Fb Yale ol dn) da cation = ie La gat ee 
. he ' aver @ Se EE lh Meek eee Std eed alae 
e028 % OW Geet be ee eR oe all be. 1 eth ahaa tillieead a ten tee ee tL 
ary ‘a Pan a a dD ad ak ae ee Te 2) Ce ede oo ie ol ted ads hee ee oe re 
BU Poh en a ’ 2 a od 2 it tor Tie Tg. 
4 ’ pO ] ° . 2 « o M er 6° ad ee ae ee he A | hs te ee bel od Ie od ee ote 2. sian 
7 A Pr) ° A a a) he Te to ie a Le Pr ed es nee odie destined amhealae Sel ena Le 
oe o : Pete @ e . e . wer Ft ed JF} -@¢-* @ wreg"? « nh a) bd cee catlenth cithetiadet atm pe to roa 
r « = 1 Oo & . wate aryee * e orgew se eo Lt at aie  iedienan ded RO ae 
« . . . te a) ‘ee ry e 5 od) eri gz: fe oh eel he« eg eee Fate al bie doe ae ten 
Pr | F ri 5 $018 ye Gage ~petroves Cad ae or Me) bP . ae phates oil a  ieotiatieoe 
y a 7 n FI A is INL yd tel ah ie eet Lette sed a te semalhaed 1! Or a eet tee Pen 
' O od a] os $e Ce i od i be ae) A i eee ee Od oF ate ae os 
ee iy t . ve ; ee « LI . oat gid "4 ous . Slt ee ee ee | a J od AE Ey el sesso 1 Ghcc deme cect a aaa 
7 ' ' ‘ A ed | i ie eS Cd ee ood bel bo De Fe I al ctte hie ae ded on ee 
7 a a Pett ee ac PTR DISET YG JCJ ek J Pee Pe , Sadi id. A ite dept! yi ate tpn rat 
' . , ae ee ‘ee: n ee Ce “¢ Ch itl eel Padi dela) ‘a ied Ieee to) ts Iek heed end e Lt 0 near os Oped on 
an) Py A ar tor Ty iN eC eee Pe e ae dd ae hel F VOe ONE TEC eOhe Mm r= g~ wee een we MaSedediee dimen. 
5 a A ee ae) ar) w A rd | od) ak ee betel beh at Rei eet ol titetietm de hee 
‘ ae Cn 0 ' o¢ hn ee 2d ual ze eo, ee he ee OT eta Dl eed 2 t pa 
' J een a A ' Cd any wr a) Pe ee Ch ood) we BS Pet MeGee ot wets? sees Ca teh tet ee onl et ne On 
f p Fi ry Py sorte oy rs ‘ rete e OE lb hd bid MIE dete Leh att Neha ot be an Inte ottarhoot hee 
. ' = ® Fy ry ed a a a Oe rae Pe el ee ot ee Od ee ee) ee ee ee ee 
‘ nd eu ea é ar] ‘ e Pl . 1¢@ a ee) sll eee fe o8 bel bedi teeiaeties Tatech Sete ates A hn Sieeeientbatimeibstteddeedt aie’ at 
O : Py r Sy ee oe ae | PT ad atl ret Sd od ST a bool he binant athe debdethan are tne Laon ty donhdeee _ 
" ' . A 6 n rs = +e n e+e pl fd Dalal ad biteieatteteh hetdeeshiel | bthteeet an afl tooo ak adikedeadaAr ion 
n Pan o a] eae o Ca a. Os i i Di Ae I Ae el, el ie nee at einen dl be olan deen a ose 
cd r ' : Te ee | ey 7 8 CO er Pn Dod Ut idea hidden penne te td Rede be olla Le Te te ee 
‘ G Sonu Je a ; ie TM btiel We hcl Ali Nite atch tirteiaaeh aentel at heed, bean a damdsdmetaah apes 
1 r fn s Py A Can) © @ ge ial od aes & te ok I betel SA ted Oh. CE a ol he te ee ve i npr ats 
Py . ‘ ' Py o ry ' Ff [) rT prey tee eR fi nt ad sitiliehthetediteha te, dot ot noe he ee 
r 8 ' 6 . A F ee Tay et te mary Pi he ee ee oO 1) a ee ee eee Oe eee ee 
U PY Pe ee re Ce 8 . cn pia foe dk de ol a of Mn Riles Be chedhcn Pein 
‘ 3 ry « a | tees a | Fa | 2 oe he) od s/ qe bel Mi deh heheh et ek tT te tag ee eee ed Ce ee 
, sg ' ¥ ap eae Le 4 ‘ i he I allah etd tie bd ch tlhe deine be Pic hein portent t 
e ‘ : e888 « ra rs A Oe bt Oh ee? ee oe | hibits adhe ten J. te chilis antes, Letetepeomennt bdees te 
Pa Pay Par) os a i « tose e i aay c Sk toed ae iat ak ltd ae phd ted ah eT Be ee 
. iD ' Fy D -o 8 Chg ® i. Ci bd , | eel ee ee ee FF Peta s@o Ly pene etatatinas id pitted 
rs © te e . oar Pa rs a ee ee ad Miedo ee et cheetah de a SN Sag ts 
e e+ ¢ ar} o en D re a | i Mid bed ted ae ee Te oe ee 
Paes, aye CN CR NC Ir aot a eli tae TRIE ot ey peetettrenaeds ancnbtenpteetn 
Py ’ ee Po Fd s re a ae | re be et w* +@e Ce Le ieee et Td Cd ed et ee Mk ee Fm fps orp Meng -ge- 
4 Py oe * Fi 5 ree eo 6 t= ee et he oe ed ad 1s * a8 gw sone ges cde cpegateademh ilies ‘aali Le ade deemed bien 
a a bat ACN ce ee Ia) Cee J a i del a ok OS SY De Lt ee ee ee aneataeth 
5 =e s Aeon AS esto ee ee ee eh Ah Od A eed htt ee Me Tey epee a ane 
5 o ' fe , Tse ee Pr ee ae ee aed A | il Sc Seok taka fe aren ict tak ect f ee note oe 
o o o , . or ar) ou ® ¢ sta gw pte O Ot lel hh tobe Se Re ete ek ee ° ©. -greura 
te id ae S F + 8 ee | a ee ee oe ee | eo: OE hdl a Od is a ee a Te iain ey eee a tate ee patndy 
. ene , a os - S 7 Se ee ceded ah ELLs seen Sel atded fi Tas el no tin dt of icbednt 
s SY ' ’ A P e oy eee ee | ad | ye ri adel oe LE et Al he el Fiat Ce ee Ae ae oe eee ee 
Pa A a ee er Toa Pr Ce) Cr ee ed a Pe eT teeter Rey Foe nT Pe why 
a oe a. ar n 1s @' oF Oe eth 8 ate Ct a = 2D ecetan ov od city © addi eee Te 
er = - A ry ence Fa) ‘ re ee ee | ee ae eb hd Oe oe to LY ee oe ae Le UI pl he Sel a 
: i A Ae = Py ry PY Pa Cr i. a de Oe ed oe ar) bid blll 4 5 deeds Reade teat bndeh the th Lids 2? onan 
‘ ’ 5 ee oe e Poy » «wt A ote Ld] Ot Le i eel ol tae oe Peat el eee 
A an 2 e864 Pore) ® 1 ar ee ee a Lo Le? Leted Se ea ee ee Le oe tt ed teat, pete 
bd er id a ® es ob A te ey ee ah ne ans ined 
0 A ¢ Cr ee ee ee 8 Prt i he i CI Ae he a ed teil hE lend ote D In hot Ann ae teh ch hated iontdat t. 
Fi eo. ; 5 er Ps ee Tet; ety § te te ptpteh ge 8 6 Fy pyegtiwicg » etn ey ee oe 
" 7 = a ee a ‘ Tee ha MC th! Ch a “4 . i | pt aba Ce ie de Bee ot Potts) eet dle ae 
n A ° "# ete fo ws e t lhe) ® ttt Y OMe VET yl vow pre 5 otet SEW Ph ae Bee 
' ae ie i SO aT NC MeO OC Re aBiats ih cath eee od ok UO Ta) es oe ee ee tO ee as neat 
D . r rr he eo Ce he i i ee Wi Ae ee ok Ke ee ey te ed 
ee . ar) eh cad i ae) ¢  —- MneFers ws Ke SAI ek San 2) nie hs Detaled eh ah ih bh 
one F . + Oe) te Pe ot be et) Pi PL Mk oth re ae a 8 pe ae np nee a ae 
. . F OC MO RD aa py Red ing) chet hs piccie ta iene a ot 
F 8 0 aoe Fo ee ie at ee ee et ep he Por res 
ees "i o o Par ek ee ey ee a ee eek ee ie tes pee) 
se » tgvete nm , a ah tl Lit ee MD ae Be te Re ha tere ot ede ella) ere gerae 
° ’ « . @ oe rn ry Pi Le oe es | Oh LP ae | io ted 44 Set Bet ho i ee a ahh ben A a sdb is dendivwcio fei one eae 
Fy . 8 et a A OD  % Le A ed ee ee ee el oe ee a ee 
: ; : , , oh ae er ee Pn ee tet tn ett te alee ey —_ 
o o P a) n ee o ee | . U A mt ae OU DM ot Bt be de ot Be ke bea he al aed C— 
q ‘ f or] o e 4 A the ie Oe I i ee ee fre ad] FP iy See ene eee el ed ete ee ee sh dl hn 
es Pea) n Pha a ee oe ° e ¢ ‘Part =e mr tok Ibori ba phn At ha tae ok) Cy e 
o oe O A were re b ts ee ee eh i ) ef Se Da a LTR Oe 
U o oan) 0 oar} ae “ a | Ce at ke LA 7%) O*Pegegrers ee ee 
« f A a@e Sd ao ah Dial ate Il hadi ok Seed tip ten cet Dood, dt bh an, Eel sh tend Lee oe 
n o re | : ! ‘Rae a oti tl kde 1a teh dhe bik | 
A ™ ae Py a ° er hee i. a ey a ae | hpaachepe tee tt te eT ee 
o | ee a to Ata ro bee ly 44 Ap ets b4 Bh nse herdh Aadond ee ann Leet net C 
as vd ae a+ a all LE od hg bn deta ahi emchan nah hdd 4m lan ie Bd doer oed ® 
4 e e ¢@ 4 ro => Sh he ee | ee tele ee te | hd 
A | Pa Pies he et ee be ee rt es tee Pe ee 
i 4 F i Oh eA dS. ke ee eee ke te eee 
p y D » n ao, Ct 1 UPS 2 PAGS MY CEPTS es CQure EeOes Pieter 
' , rs ' LU ca in | ae he | e ¥  sevyte ec ¢ led bd le bp ee dee, ks ee} i an tn li le CRA A A ha aed ok ie el 
C oan Ta aa ° ye se Sah O* meg? pak er eery ot Se he te ee oh dood tot & tho ts ee 
¥ Go. ie rr aa Sete] Bian My ty oh oc hehe 
one D ' a " - 4 eke ' Clo) a i hte Meee Je ad 
a ® ’ bs] « r i] a 5 nn 9,9 ee ae ee tad a.) 
n ay ® ee a a bell El het wt oe hE bald a) oe Thee tee | i. te) aoe 
g s & eT be ts ete hy te beets | abn) eee Lee ts 
F F ae ° a0 So adosle. nA Ser ba dy i ne Me hed healt he edo hos hea 
a ‘ ad ie) pad bed at seh is ta mie. bn Nees 
f] Fi ee A mas sh died dh oe) r OTS og Pk HHI Fe cece 
os ® | | tom » | See ek es he ee 
Py 5 e = @ bel 2 he Sih i Pie aot TD fs Pie Te fe eee 
r <6 rt or) ” om 5 hee se th ey 1 tel bah elit tie beptiectened 
Fi A 4 =. 27 ’ ih ba Dehli te teh hn Linea ba fai ethene tate 
n . aC Os 4 | i Reve es 3 Ef fete Bag ch A Nepean tctemnp 
a rs 2 a f) 7 cea De Di Po nk he mead p2 5099 SU tems 
Mas! Det ee ne bi ise. SI ae ’ 
rs ' ? _— tad, dail: a Oe bate hhh Ah ee are 
+. 8 a) r) Pig i i Lae Oo) (eo M 4 
ns r ‘ ® ee eS i a ee eS Se ee eee 


=, . = ~ 
: iain 
7 ee 
et @cC 
oC vt 
mee be) 


ci 


wm a 
oan £ 




















NAVAL POSTGRADUATE SCHOOL 


Monterey, California 





THESIS 
W 443719 


ADAPTIVE ALGORITHMS FOR TWO 
DIMENSIONAL FILTERING 


by 


Steven L. Wilstrup 


> wf hh 


september 1988 


Thesis Advisor: Murali Tummala 





Approved for public release; distribution is unlimited 


Ni 


\ 


T2424 





classified 


rity classification of this page 


REPORT DOCUMENTATION PAGE 


Report Security Classification Unclassified 1b Restrictive Markings 





Security Classification Authority 3 Distribution;Availability of Report 

Declassification Downgrading Schedule Approved for public release; distribution is unlimited. 
erforming Organization Report Number(s) > Monitoring Organization Report Number(s) 

Name of Performing Organization 6b Office Symbol 7a Name of Monitoring Organization 

val Postgraduate School (if applicable) 32 Naval Postgraduate School 


Address (city, state, and ZIP code) 7b Address (city, state, and ZIP code) 

onterey, CA 93943-5000 Monterey, CA 93943-5000 

Name of Funding/Sponsoring Organization | 8b Office Symbol 9 Procurement Instrument Identification Number 
(if applicable) 


Address (clry, state, and ZIP code) 10 Source of 10 Source of Funding Numbers i a i ss—<‘i—sO:CSCSd Numbers 


Program Element No Work Unit Accession No 


Title (include security classification) ADAPTIVE ALGORITHMS FOR TWO DIMENSIONAL Sa ENG 
Personal Author(s) Steven L. Wilstrup 


. Type of Report 13b Time Covered 14 Date of Report (year, month, day) 15 Page Count 
aster’s Thesis From To September 1988 75 


Supplementary Notation The views expressed in this thesis are those of the author and do not reflect the official policy or po- 
on of the Department of Defense or the U.S. Government. 
Cosau Codes 18 Subject Terms (continue on reverse if necessary and identify by block number) 


Id Adaptive, Filter, Wiener, LMS, RLS, Noise canceler, ALE 
Dp 
ee 


Abstract (continue on reverse if necessary and identify by block number) 

In this thesis, an adaptive two dumensional least mean squares (LMS) algonthm and a recursive least squares (RLS) al- 
rithm are developed from the one dimensional algorithms. Design of the two dimensional LMS and RLS algorithms are 
idied for accuracy based on the results of a two dimensional system identification model which was used in testing the al- 
nthms. Application of the two algonthms is demonstrated through computer simulation in which the adaptive filters are 
ployed in a noise canceler and an adaptive line enhancer and applied to an image processing problem. 





Distribution/Availability of Abstract 21 Abstract Security Classification 

unclassified ‘unlimited LJ same as report CL) DTIC users 
: Name of Responsible Individual 22b Telephone (include Area code) 22¢ Office Symbol 

all Tum a) e368 


FORM 1473,84 MAR 83 APR edition may be used until exhausted security classification of this page 
All other editions are obsolete eS 


Unclassified 


jah 


Approved for public release; distribution is unlimited. 


Adaptive Algorithms for Two 
Dimensional Filtering 


by 
Steven L. Wilstrup 
Lieutenant Commander, United States Navy 


M.B.A., University of West Florida, 1983 


Submitted in partial fulfillment of the 
requirements for the degree of 


MASTER OF SCIENCE IN ELECTRICAL ENGINEERING 
from the 


NAVAL POSTGRADUATE SCHOOL 
September 1988 


ABSTRACT 


In this thesis, an adaptive two dimensional least mean squares (LMS) algorithm and 
a recursive least squares (RLS) algorithm are developed from the one dimensional algo- 
rithms. Design of the two dimensional LMS and RLS algorithms are studied for accu- 
racy based on the results of a two dimensional system identification model which was 
used in testing the algorithms. Application of the two algorithms is demonstrated 
through computer simulation in which the adaptive filters are employed in a noise 


canceler and an adaptive line enhancer and applied to an image processing problem. 
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I. INTRODUCTION 


The area of digital signal processing has experienced a rapid growth in the last dec- 
ade. Reasons for this have been the tremendous advances in integrated circuit technol- 
ogy, and some significant developments in digital processing techniques achieved during 
this period. Included in these developments are methods which extend certain one di- 
mensional digital signal processing techniques to two dimensions. This extension is by 


no means trivial. Three significant factors must be considered: 


1. more degrees of freedom are available in a two dimensional system; this gives a 
system designer more flexibility than the one dimensional case, 


2. one dimensional problems generally involve considerably less data than two di- 
mensional problems, and 


3. the mathematical methods for handling two dimensional systems are generally less 
complete than those for one dimension. 

As the techniques for processing multidimensional data have improved, the appli- 
cations of digital signal processing have spread from one dimensional to two dimen- 
sional, to n-dimensional data. There exists many physical phenomena that inherently 
depend on two or more independent variables. In the prediction of weather and in 
seismic analysis, the data depends on more than one independent variable. Moreover, 
data originating from one dimensional processes can, in some cases, be considered two 
dimensional. For instance, data from periodic or cyclic processes can be represented as 
two dimensional arrays by using their periodicity. 

Besides the general applicability of two dimensional signal processing in the above 
cases, other areas which have experienced significant growth in recent vears include ra- 
dar, sonar, and radio astronomy. The two dimensional processing of images is also a 
very important one. Images depend on two spatial variables, and are continuous in 
nature. However, if we digitize them and assume linear models in their formation, dis- 
tortion, and recording we then have techniques that can be used in their processing. 
Depending on the applications, different processing techniques are used, notably: en- 
hancement and restoration of images, and segmentation and encoding of images. For 
details, see References 1,2,3, 4. 

This brief discussion about applications of two dimensional digital signal processing 


shows that they are found in a wide variety of fields. In this thesis, we are interested in 


extending one dimensional adaptive filtering techniques to two dimensions and then 
applying this in the area of image restoration. Once the transition from one to two di- 
mensions is understood, the extension to n-dimensional signal processing is fairly 


Straightforward. 


A. OBJECTIVES OF THE THESIS 

The use of the Wiener filter has proven to be a very powerful tool in the area of 
image restoration [Refs. 1,2, 4]. However, one disadvantage is the fact that the Wiener 
filter operates under the assumption that the image is stationary which generally is not 
the case in an image processing problem. In one dimensional signal processing the 
Wiener filter concept has provided a basis for various adaptive filtering algorithms. 
Within these adaptive algorithms, the filter possesses characteristics which can be mod- 
ified to achieve some end or objective and is usually assumed to accomplish this “adap- 
tation” automatically without the need for substantial intervention by the user. The 
adaptive filter can “learn” the signal characteristics when first turned on and thereafter 
can track changes in these characteristics. The first objective of this thesis 1s to examine 
a two dimensional Wiener model for image restoration which can then be extended to 
adaptive algorithms. Due to its relatively low computational requirements and the fact 
that it will work in a variety of signal environments, the least mean square (LMS) algo- 
rithm will be investigated first. 

The second objective of this thesis is to develop a second two dimensional adaptive 
algorithm. In this case, we will work with the recursive least square (RLS) algorithm 
which offers faster convergence than the gradient-search-type algorithms but generally 
involves a greater cost per data sample and more numerical difficulties. 

The final objective is to implement the LMS and RLS algorithms within a system 
identification model, a noise canceler, and an adaptive line enhancer. Each algorithm 
possesses several variants and various parameters which can be modified. A represen- 
tative sample of the various outputs will be examined and the results compared in order 


to see which provides the optimum solution under given conditions. 


B. ORGANIZATION OF THE THESIS 

Chapter II is designed to review the development of the one dimensional Wiener 
filter and then extend it to two dimensions where it can be incorporated into a two di- 
mensional LMS adaptive algorithm. Computer simulation of the algorithm within a 
svstem identification model is performed and the results are shown. The two dimen- 


sional RLS algorithm 1s derived in Chapter IJ] and again the results of the algorithm 


within a computer simulated system identification model are shown. Chapter IV con- 
tains the results of implementing the LMS and RLS algorithms in a noise canceler and 


an adaptive line enhancer. Conclusions concerning the results are also presented. 


Il. TWO DIMENSIONAL ADAPTIVE LEAST MEAN SQUARE 
ALGORITHM 


In this chapter, we will review the derivation of the one dimensional Wiener filter 
and then develop an algorithm to extend it to the two dimensional case. We will use the 
two dimensional Wiener filter to achieve a two dimensional least mean square (LMS) 
adaptive filter algorithm which will be demonstrated to be useful in image processing by 
applying the algorithm in a noise canceler mode and 1n an adaptive line enhancer con- 


figuration for the restoration of images corrupted by noise. 


A. ONE DIMENSIONAL WIENER FILTER 

The problem of estimating one signal from another is one of the most important in 
signal processing. In many applications, the desired signal is not available or observable 
directly. Instead, the observable signal is a degraded or distorted version of the original 
signal. The signal estimation problem is to recover, in the best way possible, the desired 
signal from its degraded replica. One typical example which we will deal with in this 
paper is an image recorded by an imaging system that has been corrupted by noise. The 
problem is to undo the noise induced distortion and restore the original image. 

This represents the classic one dimensional problem in communication theory where 
we must obtain an estimate of a signal of interest, which can be observed in the presence 
of some additive noise. In other words, the available information about the signal, 


s(m), 1s contained 1n the received signal: 
u(n) = s(n) + w(n) (Zap) 


Where w(n) is the noise. Therefore, we must process this available signal u(m) through 
an optimal processor that produces the best possible estimate of s(n). 

In order to establish a two dimensional Wiener filter, we must first develop a one 
dimensional algorithm. This task has been approached from many different directions 
[Refs. 5,6,7}]. The formulation by Haykin [Ref. 8] provides the most logical extension to 
two dimensions. First, we consider a tapped delay line filter similar to Figure 1 on page 
5. The filter consists of a set of delay elements, a corresponding set of adjustable tap 


gains or coefficients A(1), A(2),..., 447) connected to the tap inputs, and a set of adders 


u(n-2) u(n-M+2) 





Fioure 1. Tapped Delay Line Filter 


for summing the resultant outputs. The filter is driven by a random time series produc- 
ing the sequence w(7), u( — 1),..., uv — f+ 1) as the M tap inputs of the filter. 
We can express the signal produced at the filter output, y(7), by the convolution 


Sum: 


Mf 


ria) = Dau ye 1) (2.2) 


=| 


We desire a filter which in some way minimizes the difference between some desired re- 
sponse, a(m), and the corresponding value of the actual filter output. This difference can 


be denoted as 
e(m) = d(n) — y(n) (2.3) 


where e(t) is called the error signal. In Wiener theory, the filter is optimized by mini- 


nuzing the mean-square value of this error signal, e(v) . 


Let the mean-square value of the error be denoted by 
MSE = E {e*(n)} (2.4) 


where £ {.} is the expectation operator. This mean-square value is a real and positive 
scalar quantity, representing the average normalized power of the error signal, e(n). 
Substituting Equation (2.3) into (2.4) yields 


MSE = E {d’(n)} — 2 E {d(n)y(n)} + E {y’(n)}. (2.5) 


Next, substituting Equation (2.2) into Equation (2.5) and then interchanging the orders 
of summation and expectation in the last two terms, we get 


M 
MSE = E{dn)} —2 2 Ha) EA) 
M M = (2.6) 
+ YY AAU) E {un — k + un — m+ VY}. 
K= 19 n=! 


Assuming that the input signal u(n) and the desired response d(n) are jointly stationary, 
the three terms on the right-hand side of the above equation may be interpreted as fol- 


lows: 


1. The expectation E {a@(n)} is equal to the mean square value of the desired response 
d(n): 


P,= E{d’(n)}. (2.7) 


2. The expectation E {d(m)u(n — k + 1)} is equal to the cross-correlation function of 
the desired response d(n) and the input signal u(n) for the lag of k-1. We can 
therefore write the single summation term on the righthand side of Equation (2.6) 
as follows: 


M M 
DME nut ei = ) Midetk aap (2.8) 
=] — 
3. Finally, the expectation E {um—k+1)un—m+1)} 1s equal to the 
autocorrelation function of the input signal u(n) for the lag of m-k: 
r(m—k)=E {u(n—k4+ lhun—m-+1)}. (2.9) 


Accordingly, we can rewrite the double summation term on the righthand side of 
Equation (2.6) in the form 


M M 
> AAR) Elu(n — k + Du(n — m + 13 
=] m=z\ M M (2.10) 
=) > Alkyh(m)r(m — k). 
| 


Thus, substituting Equations (2.7), (2.8), and (2.10) into Equation (2.6), we find that the 


expression for the mean square error may by rewritten in the form 


M 


M M 
MSES@PR =2 > Mee aE: z KOOHO On Re (2.11) 
=] =] m=) 


By differentiating Equation (2.11) with respect to h(k) and setting it equal to zero, we 


have the following set of M simultaneous equations 


M 
p(k—1) = > holmr(m—h), k= DQM (2.12) 
m=] 


This system of M simultaneous equations 1s called the normal equations with optimum 


filter coefficients as the unknowns. With the following definitions, 
hg(1) 
ho (2) 
i= |e @ri3) 
Ay(M) 
P(9) 
pCi) 


p = (2.14) 


p(M — 1) 


r(0) r(1) r(M —1) 
r(1) r(0) xr M —2) 


R = : Bs (2.15) 


r((M—1) (M-—2) (0) 


we can rewrite the normal equations of Equation (2.12) in matrix form as 
p= Rhy. (2.16) 


This represents the discrete-time version of the Wiener-Hopf equation. 


B. TWO DIMENSIONAL WIENER FILTER 

The following derivation parallels the development by Hadhoud and Thomas, how- 
ever this research and simulation were completed separately and prior to the publication 
of Reference 9. In order to be applicable for an image processing problem, we must 
extend the formulation in the previous section to two dimensions [Refs. 10,11]. Thus is 
accomplished by developing a basic two dimensional Wiener filter as shown in 
Figure 2 on page 9. Within this filter, we use two input images, the reference array X 
and the primary input array D. The primary input image D is a two dimensional array 
which represents the ideal image plus additive noise, while the reference image X 1s noise 
that is assumed to be correlated to the noise in the primary input. Both the input arrays 
are MxM in dimension. The Wiener filter is an NxN causal FIR filter with a set of 


weights W, defined as 


W (0,0) W011)... =W(0,N—1) 
W (1,0) WAI. Wn) 
W, = | | a . (21g 


WAN —1,0) W(N-1,1) .. W(N-1,N—1) 


which minimize the mean of the squared error, e,, between the filter output and the de- 
sired input D. We designate / as the iteration number given by / = mM +n where m and 
n take on the values from 0 to M-1l. Just as in the one dimensional case, the filter out- 
put, y(77,7), 1s the convolution sum of the filter mask and the reference input X which is 


given by 


N=-1 N=1 


y(n) =) 2 W(1,k)X(m — Lyn — k). (2.18) 
i=6 =) 


During the jth iteration the input from array X is represented by X, where 


X ARRAY 





D ARRAY 


Figure 2. Two Dimensional Wiener Filter 
X(1,71) A(mn — 1) = A(m,n — N + 1) 
xj = fe: (2.19) 
X(m—N+1,2) X0m—-N+1,n-—1) 2. X(m—-N4+1n—-N+4+1) 


Using the Equations (2.17) and (2.19), Equation (2.18) can now be written as 


N-—-1 N—1 
> Wilk) X(UK) (2.20) 


y(m,n) = 
=) k=0 


for the jzh iteration. 
This output y(m,7) can now be subtracted from one pixel D(m,n) of the array D to 
produce the error signal at the jth iteration 


N-1 N-! 
= D(m,n) — 


W,(1,k)X(m — ln — k). (2.21) 


Since the purpose of the Wiener filter is to provide a set of weights which minimize 


the MSE, we can denote it as 
(2322) 


MSE = E{e?} 


where E{.} 1s the expectation operator. 
Using a mathematical derivation similar to the one dimensional case of the previous 


section, we can see that 


e? = D(mn) — >) SS Wk) Dlonyn) Xm — ln — k) 

Nay N= Ne (2.23) 
>) (Lk) WAp.g) X(m — Ln — k) X(m = pyn = q). 

f=0 k=0 p=0 g=0 


Substituting Equation (2.23) into Equation (2.22) yields 


N-} N-1 
MSE = E[D*(m,n)] — 2 W(1,k) ELD(m,n) X(m — Ln — k)] 
=0  k=0 (2.24) 


N=] N=) N=] N=! 


2 2, Y WALK) Wi(p.q) ELX(m — hn — k) Xm = pn —Q)). 
=) a= 


Defining P as the crosscorrelation matrix between the desired response D(m,n) and the 


reference input, R as the input autocorrelation matrix, and W as the optimum Wiener 


Weight matrix, we can rewrite Equation (2.24) as 


N-1 N—1 
MSE = E{D*(m,n)] - 2) 2, WoW P(Lk) 
fee) Nol N= eae 223) 
f=0 k=0 p=0 g=0 


Finally, if we minimize the MSE with respect to W,(/,k) then we have 


N=-1 N=} 
PK) =) > Wip.aR@— hq —¥) (2.26) 
p=0 q=0 


This equation is the two dimensional equivalent to Equation (2.12). The matrix form 


of Equation (2.26) 1s given as 
P= RW (2.27) 


where the elements of this equation are defined as follows 


[Ro] [Ri] LR] ... 
(A, ] [Ao] LR] .. 


R = . | ee ee es (2.28) 
CRwiid [Revg2] « TRo] 
Wo 
W, 
w=] . (2.29) 
Wy 
Po 
P, 
P=| . (2.30) 


Within the R matrix in equation (2.28) each element is a block Toeplitz matrix repres- 
ented by the equation 


where /,k,p,q range from 0 to N-I. 


C. TWO DIMENSIONAL LEAST MEAN SQUARE 

One means of obtaining an approximate solution for the optimum weight matrix, 
W, is the use of a two dimensional LMS algorithm which is depicted in Figure 3 on page 
13. This differs from Figure 2 on page 9 in that the error e, is used to update the filter 
coefficients before shifting the data window X, across the reference input for the next 
iteration. As in the one dimensional LMS algorithm, we are updating the coefficient 
matrix by adding the present matrix to a change proportional to the negative gradient 
of the error where the one dimensional instantaneous estimates of the gradient vector 
are based on sample values of the input and the error signal e(m)u(m) [Ref. 8]. For the 


two dimensional jth iteration, we define the updated matrix as 


Wis = Wi eX, (2.32) 
where 

W,., updated weight matrix 

WV, previous weight matrix 

Lu scaler multiplier controlling the rate of convergence and filter stability 


(e)(X,) estimate for the 2-D instantaneous gradient 


The previous equation can also be written as 


WV 41(k) = W(Lk) + 2 u (e) X(m— Ln — k). (2.33) 


These two equations give the two dimensional weight updating algorithm for the LMS 
adaptive filter. The algorithm we have developed here may be implemented without any 
form of matrix operations, averaging, or differentiation. 

The value of 4 may be chosen based on the desired tracking ability, steady-state 
mean square error, and convergence speed. In signal processing, there are several 
methods for determining a suitable value, however in many of these cases it requires 


knowledge of the eigenvalues and eigenvectors of the autocorrleation matrix. In image 
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Figure 3. Two Dimensional Least Nlean Square 


processing, one method which does not require these values for the computation of y 1s 
trial and error based on output image. 

An alternate method used in one dimensional design which again does not require 
a priori knowledge of the autocorrelation matrix, discussed by Bitmead and Anderson 
[Ref. 12], is the normalized LMS. Using our previous notation, this method can be ex- 
tended to two dimensions by first considering the two dimensional LMS update equation 


(2.32). In this equation, We redefine the step-size parameter, yn, for a given / and k as 


Lk) = —— 2.34 
u(i,k) p+ ohh (2.34) 


in which o, represents the normalized step size chosen between zero and two, f is an- 


other small positive constant, and o?(/,k) is one of the values from the matrix 


3; (0,0) 30,1). o(0,N—1) 
3;(1,0) (1,1). 9 of (1,N-1) 
Gj = 5 : ee ‘ (2.35) 


o;(N—1,0) of(N—1,1) ... of(N—1,N—1) 


The element required from the matrix depends upon the current values of / and & being 
used for equation (2.33). The matrix may be be intitialized with the values of 


(X(m — I,n — k))? and to update the values in a?, we use the following equation 
osi(K) = p(oj(IK)) + (1 — p)(X(m — Ln — ky’ (2.36) 


where p is a weighting factor between zero and one. This ensures that the value of yu 


does not become large enough to cause the algorithm to become unstable. 


D. IMPLEMENTATION 
1. System Identification 
In order to initially test the two dimensional least mean square algorithm, we 
established a system identification model shown in Figure 4 on page 15. Within this 
model the output of the known FIR filter, d(m,n), is defined as 


dm,n) = 4W(m,n) + .6 Wim — 1,n) 


237) 
— 3Wiman-1) + .2W(m—la-l). 


The output of the two dimensional least mean square filter, y(sm,m) , consists of a set of 
adjustable coefficients and is defined as 


y(m,n) = AO Wim,n) + Al Wom - 1,n) 


(2.38) 
+ A2W(mn—1) + A3 Wim 1jn—1). 


The adaptive filter output y(m,n) is then compared with the known system output 


d(m,n) to produce an error signal e(m,n), defined as the difference between them. 
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Figure 4. Two Dimensional System Identification 


The operation of the adaptive filter is to minimize the error signal e(r2,77) by 
providing an adaptive process for adjusting the coefficients. For this adaptive process 
we use the update equation (2.32) for the two dimensional least mean square algorithm 


developed in the previous section 
Veo = WM —peX, (Zo) 


For this model a 32x32 white gaussian noise matrix was used as the input and 
the rate of convergence can be be seen in Figure 5 on page 16 and Figure 6 on page 
17. Following 600 iterations all the coefficients had converged to within 10-3 of the ac- 
tual coefficient value and the error was 3x10-5. A computer program for this system 
identification model is given in Appendix A. 

2. Adaptive Noise Canceler 
The usual method of estimating a signal corrupted by additive noise is to pass 


the composite signal through a filter that tends to suppress the noise while leaving the 
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Figure 5. LMS System Identification-Rate of Convergence 


signal relatively unchanged. The noise canceler and adaptive line enhancer developed 
by Widrow [Ref. 13] are well-documented wavs of doing that. Adaptive noise canceling 
is a variation of optimal filtering that is highly advantageous in many applications. It 
uses an auxiliary or reference input derived from one or more sensors located at points 
in the noise field where the signal is weak or undetectable. This input is filtered and 
subtracted from a primary input containing both signal and noise. The reference input 
and the noise in the primary input are therefore correlated, and as a result the primary 


noise is attenuated or climinated by cancellation. 
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Figure 6. LMS System Identification-Rate of Convergence 


The basic noise canceling concept ts illustrated in Figure 7 on page 18. A signal 
is transmitted over a channel to a sensor that receives the signal plus noise, m. The 
combined signal and noise s+, forms the “primary input” to the canceler. A second 
sensor receives a noise m, Which is correlated in some way with the noise m. This sensor 
output provides the “reference input” to the canceler. The noise », is filtered to produce 
an output y that is a close replica of m,. This output is subtracted from the primary 


put s+ to produce the system output s +n, — y. 
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Figure 7. Adaptive Noise Canceler 


In the system shown in Figure 7, the reference input is processed by our two 
dimensional least mean square filter. Thus the filter operates under changing conditions 
and will readjust itself continuously to minimize the error signal. 

The computer program which simulates this noise canceling model is provided 
in Appendix B. We will discuss the results and conclusions concerning various simu- 
lations in Chapter 4. 

3. Adaptive Line Enhancer 

A special case of adaptive noise canceling ts when there is only one signal x(v) 
available which is contaminated by noise. In such a case, the signal x() provides its 
own reference signal j}(), which is taken to be a_ delayed replica of 
x(1): y(n) = x(n — A), as shown in Figure 8 on page 19. The adaptive filter will respond 
by canceling any components of the main signal x(/) that are in any way correlated with 
the secondary signal y(77) = x(m— A). Suppose the signal x() consists of two compo- 
nents: a narrowband component that has long-range correlations and a broadband 
component which will tend to have short-range correlations. One of these could repre- 


sent the desired signal and the other an undesired interfering noise. Suppose that the 
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Figure 8. Adaptive Line Enhancer 


delay A ts selected so that it falls between the correlation lengths. Since A 1s longer than 
the eflective correlation length of the broadband component, the delayed replica will be 
entirely uncorrelated with the broadband part of the main signal. The adaptive filter 
will not be able to respond to this component. On the other hand, since A 1s shorter 
than the correlation length of the narrowband component, the delayed replica that ap- 
pears in the secondary input will be correlated with narrowband part of the main signal 
and the filter will respond to cancel it. Note that 1f A is selected to be longer than both 
correlation lengths, the secondary input will become uncorrelated with the primary input 
and the adaptive filter will turn itself off. In the opposite case, when the delay A ts se- 
lected to be less than both correlation lengths, then both components of the secondary 
signal will be correlated with the primary signal, and therefore the adaptive filter will 
respond to cancel the primary x() completely. The computational algorithm for the 


adaptive line enhancer 1s shown in the following three equations: 


M M 
x(n) = > Aor) (n —m)= > Pl tt)2x( —m-— A) (2.40) 


m= m=) 


e(n) = x(n) — x(n) (2.41) 
h(n +l =h,(n) + 2pe(n)x(n—m—A) m=0,1,2,....M (2.42) 


For this model we also developed a computer simulation which incorporates our 
two dimensional LMS algorithm and it 1s provided in Appendix C. The results and 


conclusions will again be discussed in Chapter 4. 
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IW. TWO DIMENSIONAL ADAPTIVE RECURSIVE LEAST SQUARES 
ALGORITHM 


In the previous chapter, under the LMS algorithm, the available data samples were 
used in order to attempt to move the current estimate of the impulse response to the 
optimum value. This approach has the advantage of being simple to implement but 
carries with it the disadvantages that it can be slow to approach the optimal weight 
vector and, once close to it, will usually fluctuate around the optimal vector rather than 
actually converge to it due to the effects of approximations made in the estimate of the 
performance function gradient. 

To overcome these difficulties, we examine another approach in this chapter which 
uses the input data in such a way as to ensure optimality at each step. This alternative 
algorithm is based on the exact minimization of the least square criteria. This algorithm 


is known as recursive least squares (RLS). 


A. ONE DIMENSIONAL RECURSIVE LEAST SQUARES 

As in the LMS case, the one dimensional RLS algorithms is developed in several 
different ways [Refs. 5, 8 ]. Orfanidis [Ref. 6] provides a derivation which we will con- 
sider prior to our extension to two dimensions. The tapped delay line shown in 
Figure 9 on page 22 will provide the reference for the following discussion. We begin 


by replacing the LMS estimation criteria of MSE = E[e?] by 


MSE = » eh k = 0,1,...,n (3.2) 
K=0 
where 
e(k) = x(k) — x(k) (3.3) 


and x(k) is the estimate of x(k) produced by the Mth-order Wiener filter 


M 
£ = ) hm) (k—m) (3.4) 


m=0 


za 
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Figure 9. One Dimensional RLS Reference (TDL) 


Substituting equation (3.2) into (3.1) and setting the derivative with respect to h to zero 


we find the least-squares analogs of the orthogonality equations 


to 
to 


OMSE . 
SE = - 2) ell y(k) = 0 (3.5) 


which may be rewritten in their normal equation form as follows 


> — h’y(k))y(k) = 0 (3.6) 
k=0 
> weve | h = ) xy(W (3.7) 
k=0 =( 
Defining the quantities 
A Hk) yk) 3.8 
(n) Dv) ¥W (3.8) 
n) = k) v(k 3.9 
r() D2) (3.9) 


we can then write the normal equation as R(m)h = r(n), with solution h = R(x)"'r(n). 


Note that the n-dependence of R(m) and r(m) depend on the current time 1, therefore, 
h(n) = R(n)~'r(n) = P(n)r(n) (3.10) 


where P(z) = R(n)-'. These are the least squares analogs of the ordinary Wiener sol- 
ution, with R(m) and r(m) playing the role of the covariance matrix E[y(m)y(n)"] and 
cross-correlation vector E[.x(n)y(n)], respectively. The RLS algorithm is obtained by 
writing equation (3.10) recursively in m and then using the following matrix inversion 
lemma [Ref 5] 


Cen = Aro Ae BDA B+ Cy) DAS (3.11) 
we get the update equation for the P matrix 


P(n — l)y(n)y(n)‘ P(n — 1) | 


Pin) = P(n—-1) - 
7 ie 1 + y(n)" P(n— 1)y(n) 


(3.12) 
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Using the quantities in equations (3.8) and (3.9) to satisfy the recursions yields 
Rin) = R(n—1) + y(n)y(n)" (3.13) 
rim) = r(n—1) + x(n)y(n) (3.14) 


and substituting equations (3.12) and (3.14) into (3.10) after some mathematical manip- 


ulations we find 
h(z) = h(n —1) + P(r) y(n) e(n) (3.15) 


which differs from the LMS alogorithm by the presence of the P(n), vice p, in front of 
the weight correction term. Since P(m) = R(n)-' 1s a measure of the covariance matrix 
ELy(n)y(n)7J, the presence of R(m)-! makes the RLS algorithm behave like Newton's 


method, and hence has fast convergence properties. 


B. TWO DIMENSIONAL RECURSIVE LEAST SQUARES 

Although it 1s discussed briefly by Wellstead and Caldas Pintos [Ref 14], limited 
information 1s available in the open literature in this area. In order to develop a two 
dimensional recursive least square (RLS) algorithm we will extend the one dimensional 
adaptive algorithm developed in the previous section in a method similar to that used for 
the LMS algorithm. 

As in Chapter 2, using Figure 10 on page 25 as a reference we see that the basic 
filter has two input images. The primary two dimensional input array, D, is the ideal 
image plus additive noise. The reference image X is the noise array. Each array is of 
dimension M by M. The filter mask, W, 1s N by N. 

The one dimensional RLS algorithm, equation (3.15), is the same as the one di- 
mensional LMS algorithm with the exception that P(m) replaces yw. Therefore, if a 
method of developing a P matrix for the two dimensional case can be devised we can 
then use an algorithm similar to the two dimensional LMS algorithm to update the filter 


coefficients. First, let 
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Figure 10. Two Dimensional RLS Reference 


to 
csr 


P,(0,0) P(0,1) 7 P(0,N* 1) 
P, = . 4 . (3.16) 
P(N? —1,0) P({N*—1,1) ~ P(N* -1,N? -1) 


which represents the P matrix on the jth iteration where / is defined as j = mM + n. 


As in the two dimensional LMS algorithm, we let the filter coefficient matrix be given 


as 


W (0,0) W011)... W0,N-1) 
W,(1,0) Wj). = WN = 1) 


W, = | | i | aan 
W(N—1,0) WAN—1,1) ... WAN-1,N—1) 


and the input data window as 
A(m,n) X(myn — 1) se X(m,n — N+ 1) 
iT. | . | Gu) 
Xm-N4 1,n) Mm-N+ In-1) Xn —N-+1n—N+ 1) 


In order to establish an update equation for the P matrix in which each element will 
have the proper dimensions, we use equations (3.17) and (3.18) and make the following 


transformations; the filter mask, W,, equation (3.17) is transformed into a vector defined 


as 
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W,(0,0) 
(0,1) 


W(0,N — 1) 
W(1,0) 
W141) 


WW, = (3.19) 


W(1,N — 1) 


WN — 1,0) 
WN — 1,1) 


W(N —1,N—1) 


and we transform the X;, matrix, equation (3.18) into 


X(m,n) 
X(m,n — 1) 


X(myn — N+ 1) 
WX, = . (3.20) 


X(m — N + 1,n) 
A(m— N+ 1,n—1) 


X(m—- N+ 1n—-N+4+1) 
Using the quanitities defined above, we can then write the P matrix update equation as 


ps 


P._, CO ean, 


P Te 
1+ (XX/)P)_(XX) 


As discussed earlier, the one dimensional RLS weight updating algorithm differs 
from the LMS in that it contains the P matrix while the LMS algorithm contains yu. This 
must be considered in two dimensions since P is a matrix and yw is a scalar. Therefore 


using the previous equations, we can obtain the update to the filter coefficients via 
WW. = WW, — (P)e(XX) (3.22) 
where e, is given as 


e, = D(m,n) — ylm,n) (3.23) 


y(m,n) is the convolution sum of the image in the reference input with the filter window 


N—1 


N-] 
aan > 2, WAUk) X(1,K) (3.24) 
==() c==() 


and the value P, is updated by equation (3.21). 
The cost function of the RLS criterion could be modified to include a windowing 


of the input data as follows 


MSE = da e*(k) (3.25) 
k=0 


This modification results in the following modified P matrix update: 


| Ps ( XOQ OG as 
P= | P= = = (3.26) 
1 + (XX; )P,_y(XX)) 


where q, the averaging or “forgetting factor” is a positive constant. It is usually chosen 
to be slightly less than one, thereby diminishing the contribution of the “older” data. 
The problem of data nonstationarity is the main reason for introducing this type of 
weighting factor.[Ref. 6] It should be noted that the usual RLS algorithm is attained 


when g = | and that no changes in the amount of computation 1s required. 
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C. IMPLEMENTATION 
1. System Identification 
Following the mathematical development of the two dimensional recursive least 
squares algorithm, we implemented and tested it as was previously done with the two 
dimensional least mean square algorithm. Our initial testing was done in a system 
identification model as shown in Figure I}. As in the two dimensional LMS case, we 


used a known filter with the following output equation 


d(m,n) = .4V(m,n) + .61V¥ On — I,n) 


— 31 (niyn— 1) + .20KOn— 1,n-1) (ey) 
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Figure 11. Two Dimensional RLS System Identification 
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The two dimensional adaptve recursive least square filter output was given as 


y(m,n) = AO W(m,n) + Al W(m— 1,n) ; 
+ A2W(mn—1) + A3 Wim — In 1) (3.28) 


The error signal was generated by taking the difference between the desired (known) 
signal, d(m,n), and the filter output, y(m,”). By using equations (3.21) and (3.22) to up- 
date the P matrix and the filter coefficients, we caused the filter weights to move toward 
their optimum values and the error signal to approach zero. 

The computer program for this model is provided is Appendix D. For this sys- 
tem, a 32x32 white gaussian noise matrix was generated for the input to both the known 
filter and the RLS adaptive filter. The rate of convergence is shown in Figure 12 on 
page 31 and Figure 13 on page 32 and as can be seen after approximately 70 iterations 
each of the coefficients had converged to with 10-3 of the actual value and the error was 
Ge 

2. Adaptive Noise Canceler/Adaptive Line Enhancer 

These two svstems were both discussed in Chapter 2 with regard to implemen- 
tation of the two dimensional LMS algorithm. In order to further test the RLS algo- 
rithm we also implemented it within the noise canceler and the adaptive line enhancer. 
The computer programs which were used for the simulation are given Appendix E and 
Appendix F, respectively. The simulation results and discussion will be provided in 


Chapter 4. 
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Figure 12. RLS System Identification Rate of Convergence 
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Figure 13. RLS System Identification Rate of Convergence 


IV. RESULTS AND CONCLUSIONS 


In this chapter, we present the experimental results from Gur implementation of the 
wwe dimensional LMS algorithm derived in Chapter H] and the two dimensional RLS 


algorithm derived in Chapter Lil. The two algorithms were used in a noise canceler and 


i 


an adaptive line enhancer, as discussed in Chapter I], to solve an image restoration 


Figure [4 shows the original image of a house which is cornprised of 128 by 128 
eight bit pixels. The image has a mean value of 131.68 and a variance of 3194.4. In 


x 
Caer 


Figure 15 on page 34, we have corrupted the original image with additive white gaussian 
{ 


zero mean and variance 1600). 








RTL FE PEF FO Eee TEES re TEEPE Ed Pal hE Pal ed Pal OIE Hh tel OP PIS OOTP IE Pl Pol TEED POS aI Be Pl Ps Pot Dal PO PE PN I a al I BT OL RE OT PE ek Bal ate Al 





Figure 14. Original Image 
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Figure 15. Image plus Additive Noise 


A. NOISE CANCELER RESULTS 

In order to solve our image restoration oroblem, we initially implemented toth the 
LMS and the RLS algorithms in a noise canceler using a two by two filter matrix and 
a value of 70x10° for uw in the LMS algorithm. 

Figure 16 on page 35 shows the output when the LMS algorithm 1s processed 


through the image for one pass: this results in 10,384 updates of the fliter corresponding 


c 


mys 


to the 16,384 pmels. Figure 17 on page 36 is the results following two passes throug 
the image. No significant improvernent was reakzed with more than two passes. 

implementing two by two RLS adaptive filter in the noise canceler produced very 
favorable results after one pass through the image. These results are shown in 
Figure 18 on page 36 and compare weil with the criginal image and the LMS output 
after tWo passes. 

For the LMS algorithm our best results were achieved with a value of 70x10-° for 
u, however in order to demonstrate the effect of changing this value Figure 19 on page 


37 represents the output with « = 34x10° and Figure 20 on page 37 is produced Gy a 


spatially varving normalized 4 as discussed in Chapter [I]. I: can be seen that VAarlous 


w2 
an 


values of u, different features within the image were restored at different levels. In Fig- 
ure 16 on page 35, the edges and more detailed segments of the image were restored 
better than in Figure 19 on page 37, however areas of simular contrast were not as well 
restored. When using the normalized uw the mean value more closelv aprroached the 
original mean, however the variance was reduced. 

Finally, imereasing the number of coefficients, using a three bv three filter matrix in 
the LMS algorithm we obtained the results shown in Figure 21 on page 38 in one pass. 
This filter showed no significant improvement m the variance compared ta the two bv 
two filter; however, 1 did improve with respect to the mean. 


Table | on page 38 shows the restoration results comparing the mean and variance 


of the various cuiputs with that of the original image and the image plus noise. 





Figure 16. Restored Image: Noise Canceler/LMS (One pass) 





Figure 17. Kestored Image: Noise Canceler/LMS (Two passes} 











i ne a 
FETIP EELE ne 
Py 






Aah P at EDL he fo tn Mate md hy nes 
Arh 





fire Botner 


Bo Pr DSP PEI OLR NER oy 3: 
ee LEIRPY IELTS 





fe " 


SCR TIRBRPR Deo Ge fstabs 


‘ , oy 









Adin timetan fa tin Lepetie tion che etna tin ape D ln fae tien Ai han tar in Lane hine fee en foe thn neon cen ite Dm Bhar tin tiny So, Ue re Teen tant ane an te ae hen faenchen ti ne cen Tcor ete De Pha hn ch fe Lon mectanetin enn enn ee enfin dan hy ete ean Lenten daeetin ane tan cee ane lon enced ee cnn ene inlen Anchen re ca endiensieeieeeien elie ed 





‘joure 18. Restored image: Noise Canceler/RLS (One pass 





Figure 20. Restored Image: Noise Canceler/LMS (normalized mu) 





Figure 21. Restored Image: Noise Canceler/LMS (3 by 3) 


Table 1. NOISE CANCELER RESULTS 


One pass LM S(different 120.12 1940.4 
mu) 

One pass 126.53 1500.7 
LMS(normalized mu) 


One pass 3x3 LMS iS 2017.3 
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B. ADAPTIVE LINE ENHANCER RESULTS 

We next implemented our two algorithms in an adaptive line enhancer. In order to 
obtain a reference input frorn the primary input, we used a two dimensional delay oper- 
ator of (m-i.nj. As in the noise canceler we again used a two bv two filter matrix and 
this produced the results shown in Figure 22 on page 39 and Figure 23 on page 40 for 
the one pass LMS and the one pass RLS, respectively. No significant improvement was 
noted in the two pass LMS compared to the one pass LMS, 

Tabie 2 on page 40 shows the comparison of cach output mean and vanance with 


the mean and variance of the original image and the image plus noise. 
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Figure 22. Restored Image: Adaptive Line Enhancer/LMS (One pass) 
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Figure 23. Restored Image: Adaptive Line Enhancer/RULS (One pass} 


Table 2. ADAPTIVE LINE ENHANCER RESULTS | 
| TYPE OF FILTER | MEAN sf VARIANCE 
Original fina age i31.68 3194.4 
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ige plus Noise $34.25 ri6a77 


One pass 2xZ L. MS Bees L 2531, 8 
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C. CONCLUSIONS 


in this thesic we have extended a one dirnensional LMS and 4@ one dimensional KLS 


adaptive tes onthm to two dimensions. As we examine the results we see that many of 


ihe one dimensional comparisons between LMS and RLS are applicable to two dimen- 


SiQins: 


ie 


the rate of convergence of the RLS algorithm is in general faster than that of the 
LMS algornthm by an order of magnitude 


this s oo performance of the RLS algorithm is atrained at the expense of a large 
increase in the computational complexity 


3. there are no approximations within the derivation of the RLS algorithm, therefore 
as the number of iterations approach infinity, the least squares estimate approaches 
the optimum Wiener value 


4. in the RLS algorithm the correction which is applied is computed using all past 
data whereas the LMS algorithm uses only the instantaneous sample and the error 
signal; this is not necessarily an advantage in image processing. 

The objectives of this thesis were successfully completed. Some suggestions for fu- 
ture work include: (i) to examine these two dimensional algorithms in other applica- 
tions, applying them to areas other than image processing, (11) to extend the two 
dimensional LMS and RLS algorithms to n-dimensions and implement them, (iii) to 
analyze the extension of these one dimensional algorithms to two dimensions in the fre- 
quency domain. Multidimensional digital signal processing is being applied in many 


areas today, however the potential for future growth and applications appears unlimited. 
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APPENDIX A. LMS SYSTEM IDENTIFICATION PROGRAM 


SYSTEM IDENTIFICATION OF AN ADAPTIVE TWO DIMENSIONAL 
LMS FILTER VS. A KNOWN TWO DIMENSIONAL FILTER 


REAL MU 
DIMENSION AO(0: 32,0: 32) ,A1(0: 32,0: 32) ,A2(0: 32,0: 32), 


* A3(0: 32,0: 32),E(0: 32,0: 32),D(0: 32,0: 32) ,W(-1: 32,-1: 32), 
* Y(0: 32,0: 32) 


dekreters VARIABLE DEFINITIONS #2 %veve% 
A=NOISE AMPLITUDE 
S=NOISE VARIANCE 
AM=NOISE MEAN 
IX=SEED 
AO,A1,A2,A3=FILTER COEFFICIENTS 
MU=SCALING FACTOR WITHIN THE LMS ALGORITHM 
W=WHITE GAUSSIAN NOISE GENERATED BY SUBROUTINE 
Y=FILTER OUTPUT 
=DESIRED OUTPUT 
E=ERROR(D - Y) 


tedekeacece INITIAL VALUESG* xe 
A=1. 

MU=. 01 

A0(0,0)=0. 

A1(0,0)=0. 

A2(0,0)=0. 

A3(0,0)=0. 

S=1. 

AM=0. 

IX=65539 


weet QPEN FILE AND FILL INPUT MATRIX WITH WHITE 
GAUSSIAN NOISE este veve te deve 
OPEN (UNIT=80,FILE='SYSID2 DATA' ,STATUS='NEW' ) 
DO 10 M=-1,31 
DO 20 N=-1,31 
IF((M. LT. 0). OR. (N. LT. 0)) THEN 
W(M,N)=0. 
ELSE 
CALL GAUSS(IX,S,AM,V) 
WCM,N)=V*¥A 
ENDIF 
CONTINUE 
CONTINUE 


#eeCOMPUTE THE UPDATED FILTER COEFFICIENTS 
USING THE TWO DIMENSIONAL LMS ALGORITHM**%* 
DO 40 K=0,31 
DO 50 J=0,31 
D(K,J)=. 4*W(K,J)+. 6*W(K-1,J)~-. 3*W(K,J-1)+. 2*W(K-1,J-1) 
Y(K,J)=A0(K,J)*W(K,J)+A1(K,J)*W(K-1,J)+A2(K,J)*W(K,J-1)+ 
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50 


aa 
Oo 


15 


70 
60 


50 


A3(K,J)*W(K-1,J-1) 

E(K,J)=D(K,J)-Y(K,J) 
AQ(K, J+1)=A0(K,J)+MU*E(K,J)*W(K,J) 
A1(K,J+1)=A1(K,J)+MU*E(K,J)*W(K-1,J) 
A2(K, J+1)=A2(K,J)+MU*E(K,J)*W(K,J-1) 
A3(K, J+1)=A3(K, J)+MU*E(K, J)*W(K-1,J-1) 

CONTINUE 

AO(K+1,0)=A0(K,J) 

Al(K+1,0)=A1(K,J) 

A2(K+1,0)=A2(K,J) 

A3(K+1,0)=A3(K,J) 

CONTINUE 


Heke KOUTPUT RESULT Sitters 
L=0 
DO 60 K=0,31 
DO 70 J=0,31 
L=L+1 
PRINT 15, L,E(K,J),A0(K,J),A1(K,J) ,A2(K,J) ,A3(K,J) 
WRITE (80,15) L,E(K,J),A0(K,J),A1(K,J),A2(K,J),A3(K,J) 
FORMAT (' ',1X,14,2X,F8.5,2X,F8.5,2X,F8.5,2X,F8.5, 
2X,F8. 5) 
CONTINUE 
CONTINUE 
STOP 
END 


SUBROUTINE GAUSS(IX,S,AM,V) 
A=0. 0 

DO 50 I=1,12 

CALL RANDUCIX,IY,Y) 

IX=1Y 

A=A+Y 

V=(A-6. 0)*S+AM 

RETURN 

END 


SUBROUTINE RANDUCIX,IY,YFL) 
TY=IX*65539 

Bimey )5,6,6 
TY=ITY+2147483647+1 

YFL=IY 

YFL=YFL*. 4656613E-9 

RETURN 

END 
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APPENDIX B. LMS ALGORITHM IMPLEMENTED IN A NOISE 


CANCELER 


THIS IS A VAX/VMS FORTRAN PROGRAM THAT IMPLEMENTS A TWO 
DIMENSIONAL 2X2 ADAPTIVE LMS FILTER WITHIN A NOISE 
CANCELER 


INTEGER M,N,K,J 
BYTE BC0: 127) ,BYTEECO?: 127) {BY TEUC0: t2gee 
INTEGER*4 INTEC0: 127 ,0: 127) , INTUCO: 127,0: 127), 
LEU 
REAL*4 MU,A,FMINE, FMAXE , FMINU, FMAXU 
REAL*4 AAC 0: 3),EC0: 127,0: 127) ,UC0: 127 ,0: 127) 
RBAL*4 WC-1:127,-1: 127), YC0: 127,03 127) oi 127 0 


Wieder VARTABLE DEFINITIONS rrr 
A=NOISE AMPLITUDE 
S=NOISE VARIANCE 
AM=NOISE MEAN 
AA=FILTER COEFFICIENTS 
IX=SEED | 
MU=SCALING FACTOR WITHIN THE LMS ALGORITHM 
IM=INPUT IMAGE 
W=WHITE GAUSSIAN NOISE GENERATED BY SUBROUTINE 
U=SIGNAL PLUS NOISECIM + W) 
Y=FILTER OUTPUT 
E=ERROR(U - Y) 
FMINE , FMAXE , FMINU , FAMXU=PARAMETERS 
TO BE USED TO CONVERT DECIMAL DATA TO BYTE DATA 


Fevedkeskeseve de T NITIAL VALUE GO rerededeveuedeve ve 
A=1. 

MU=35. E-9 
AA(0)=0. 
AA(1)=0. 
AA(2)=0. 
AA(3)=0. 

S=40, 

AM=0. 
IX=65539 
FMINE=1. E+10 
FMAXE=-1. E+10 
FMINU=1. E+10 
FMAXU=- 1. E+10 


wierererrQPEN AN IMAGE FILE, CONVERT THE BYTE DATA INTO INTEGERS 
AND THEN PLACE THESE VALUES IN A MATRIX*#X?eivevex 
OPEN (UNIT=1, NAME ='HOUS1G.DAT', TYPE ='OLD', ACCESS = 
‘DIRECT’, RECORDSIZE=32, MAXREC=128) 
DO 100 M=0,127 
READ(1'M+1) (B(J),J=0,127) 
DO 110 N=0.127 


+r 


P10 
100 


50 
40 
C 
C 
C 


550 
500 


IF(B(N). LT. 0) THEN 
IM(M,N)=B(N)+256 
ELSE 
IM(M,N)=B(N) 
ENDIF 
CONTINUE 
CONTINUE 


HietereIeAKADD WHITE GAUSSIAN NOISE TO THE IMAGE AND SET THE 
VALUES OUTSIDE THE IMAGE TO ZERO**%*xr7¥% 
DO 10 M=-1,127 
DO 20 N=-1,127 
eeCheign 0). OR. CN. LT. 0) ) THEN 


W(M,N)=0. 
IM(M,N)=0. 
ELSE 
CALL GAUSS(IX,S,AM,V) 
WOM,N)=V*A 
ENDIF 
U(M,N)=IM(M,N)+W(M,N) 
CONTINUE 


CONTINUE 


sek USE THE LMS ADAPTIVE ALGORITHM TO UPDATE 
THE FILTER COEFICIENTS*6iiiesesriek 
DO 40 K=0,127 
DO 50 J=0,127 
Y(K,J)=AA(0)*W(K,J)+AA(1)*W(K-1,J)+ 
AA(2)*W(K,J-1)+AA(3)*W(K-1,J-1) 
E(K,J)=U(K,J)-Y(K,J) 
AA(0)=AA(0)+MU*E(K,J)*W(K,J) 
AA( 1)=AA( 1)+MU*E(K,J)*W(K-1,J) 
AA( 2)=AA( 2)+MU*E(K,J)*W(K,J-1) 
AA(3)=AA(3)+MU*E(K,J)*W(K-1,J-1) 
CONTINUE 
CONTINUE 


week CHANGE SIGNAL PLUS NOISE AND ERROR 
OUTPUT INTO BYTE DATA AND THEN WRITE 
LO) on FILExt vevedveve 

OPEN (UNIT=2 ,NAME=' ERROR. DAT' , TYPE='NEW' , ACCESS= 


* ‘DIRECT’ ,RECORDSIZE=32 , MAXREC=128) 


OPEN (UNIT=3 ,NAME=' SIGNOISE. DAT' , TYPE='NEW' , ACCESS= 


ae ‘DIRECT’ , RECORDSI ZE=32 , MAXREC=128) 


DO 500 I=0,127 
DO 550 J=0,127 
IF(E(1I,J). LT. FMINE)THEN 
FMINE=E(I,J) 
ENDIF 
IF(E(I,J). GT. FMAXE) THEN 
FMAXE=E(I,J) 
ENDIF 
CONTINUE 
CONTINUE 
IF (FMINE. LT. 0) THEN 
FMAXE=FMAXE -FMINE 
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650 


600 


750 


700 


450 


400 


ENDIF 
DO G0Cn1—0 9 lo 
DO) 650 J—omler 
IF(FMINE. LT. 0)THEN 
E(I,J)=E(1,J)-FMINE 
E(I,J)=E(1,J)*255. /FMAXE 
ELSE 
E(I,J)=(E(I,J)=FMINE)*255. /FMAXE 
ENDIF 
IE=NINT(E(1I,J)) 
IF(IE.GT.127) THEN 
IE=IE-256 
ENDIF 
BYTEE(J)-IE 
CONTINUE 
WRITE (2'I+1) (BYTEE(N) ,N=0,127) 
CONTINUE 
DO 700 I=0,127 
DON 7500 J—Ommer, 
IF(U( I,J). LT. FMINU) THEN 
FMINU=U(I ,J) 
ENDIF 
IF(UC I,J). GT. FMAXU) THEN 
FMAXU=U(I,J) 
ENDIF 
CONTINUE 
CONTINUE 
IF (FMINU. LT. 0) THEN 
FMAXU=FMAXU-FMINU 
ENDIF 
DO" 400s Cnni2y, 
DO 450) I-07 
IF(FMINU. LT. 0)THEN 
UCI,J)=U(1,J)-FMINU 
UCI, J)=U(1,J)*255. /FMAXU 
ELSE 
UCI, J)=(UCI,J)=FMINU)*255. /FMAXU 
ENDIF 
IU=NINT(UCI,J)) 
IF(IU. GT. 127) THEN 
IU=IU-256 
ENDIF 
BYTEU(J)-IU 
CONTINUE 
WRITE (3'I+1) (BYTEU(N) ,N=0,127) 
CONTINUE 
CLOSE (UNIT=2) 
CLOSE (UNIT=3) 
STOP 
END 


SUBROUTINE GAUSS(IX,S,AM,V) 
AMP=0. 0 
DO 50 I=1,12 

RANDOM=RAN( IX) 
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50 


AMP=AMP+RANDOM 
CONTINUE 
V=( AMP-6. 0)*S+AM 
RETURN 
END 
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LINE ENHANCER 


THIS IS A VAX/VMS FORTRAN PROGRAM THAT IMPLEMENTS A TWO 
DIMENSIONAL 2X2 ADAPTIVE LMS FILTER WITHIN AN ADAPTIVE 
LINE ENHANCER 


INTEGER M,N,K,J 

BYTE B(0: 127) ,BYTEY(0: 127), BYTEU(0: 127) 

INTEGER*4 INTY(0: 127,0: 127), INTU(0: 127,0:127),IY,IU 
REAL*4 MU,A,FMINY,FMAXY , FMINU, FMAXU 

REAL*4 AA(0: 3),E(0: 127,0:127),U(0: 127,0: 127) 

REAL*4 W(-1:127,-1:127),Y(O: 127, 0: 127), IM(0: 127,0: 127) 


WwiekkeKVARTABLE DEFINITIONS ¥%X v2r% 
A=NOISE AMPLITUDE 
S=NOISE VARIANCE 
AM=NOISE MEAN 
AA=FILTER COEFFICIENTS 
IX=SEED 
MU=SCALING FACTOR WITHIN THE LMS ALGORITHM 
IM=INPUT IMAGE 
W=WHITE GAUSSIAN NOISE GENERATED BY SUBROUTINE 
U=SIGNAL PLUS NOISECIM + W) 
Y=FILTER OUTPUT 
WW=DELAYED VERSION OF SIGNAL PLUS NOISE(CU) 
E=ERROR(U - Y) 
FMINY , FMAXY , FMINU , FAMXU=PARAMETERS 
TO BE USED TO CONVERT DECIMAL DATA TO BYTE DATA 


wk INITIAL VALUES * kkk % 
A=1. 
MU=35.E-9 
AA(0)=0. 
AA(1)=0. 
AA(2)=0. 
AA(3)=0. 

S=40. 

AM=0. 
TX=65539 
FMINU=1. E+10 
FMAXU=-1. E+10 
FMINY=1. E+10 
FMAXY=-1. E+10 


wkieteeedeekOPEN AN IMAGE FILE, CONVERT THE BYTE DATA INTO INTEGERS 
AND THEN PLACE THESE VALUES IN A MATRIX %e%iro% 
OPEN (UNIT=1, NAME ='HOUS1G. DAT', TYPE ='OLD', ACCESS = 
DIRECT’, RECORDSIZE=32, MAXREC=128) 
BO 100 M=0,127 
READ(1'M+1) (B(J),J=0,127) 
BO melOmN=0027 
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IFC BCN). LT. 0)THEN 
IM(M,N)=BC(N)+256 


ELSE 
IM(M,N)=B(N) 
ENDIF 
ae CONTINUE 
100 CONTINUE 
C 
C Fee ADD WHITE GAUSSIAN NOISE TO THE IMAGE AND SET THE 
C VALUES OUTSIDE THE IMAGE TO ZERO?eideddcdss 
DO 10 M=-1,127 
DO 20 N=-1,127 
IF ((M. LT. 0). OR. (N. LT. 0))THEN 
W(M,N)=0. 
IM(M,N)=0. 
ELSE 
CALL GAUSS( 1X,S,AM,V) 
W(M,N)=V*A 
ENDIF 
U(M,N)=IM(M,N)+W(M,N) 
WW(M,N)=U(M-1,N) 
20 CONTINUE 
10 CONTINUE 
E 
C ‘eeteeIUSE THE LMS ADAPTIVE ALGORITHM TO UPDATE 
e THE FILTER COEFICIENTS *ietvedevedescoex 


DO 40 K=0,127- 
DO 50 J=0,127 
Y(K,J)=AA(0)*WW(K,J)+AA(1)*WW(K-1,J)+ 
os AAC 2)*WW(K,J-1)+AA(3)*WW(K-1,J-1) 
E(K,J)=U(K,J)-Y(K,J) 
AA(0)=AA( 0 )+MU*E(K,J)*WW(K, J) 
AA(1)=AA(1)+MU*E(K,J)*WW(K-1,J) 
AA(2)=AA( 2)+MU*E(K,J)*WW(K,J-1) 
AA(3)=AA(3)+MU*E(K,J)*WW(K-1,J-1) 
50 CONTINUE 
40 CONTINUE 
C 


c weFeeAe eK CHANGE SIGNAL PLUS NOISE AND FILTER OUTPUT 
G INTO BYTE DATA AND THEN WRITE TO A FILE*tevervevs 
OPEN (UNIT=3,NAME='SIGNOISE. DAT’ , TYPE='NEW' , ACCESS= 
% 'DIRECT’ , RECORDSIZE=32 , MAXREC=128) 
OPEN (UNIT=4,NAME='FILTERED. DAT’ , TYPE='NEW' , ACCESS= 
¥ "DIRECT' , RECORDSIZE=32 , MAXREC=128) 
DO 700 I=0,127 
DO 750 J=0,127 
IF(UCI,J). LT. FAINU)THEN 
FMINU=U( I,J) 
ENDIF 
IF(U(1,J). GT. FMAXU)THEN 
FMAXU=U( I,J) 
ENDIF 
750 CONTINUE 
700 CONTINUE 
IF (FMINU. LT. 0) THEN 
FMAXU=FMAXU-FMINU 
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450 


400 


850 


800 


950 


900 


ENDIF 
DO 400 I=0,127 
DOr4s One, 
IF(FMINU. LT. 0)THEN 
UCI, J)=U(1,J) -FMINU 
UCI,J)=U(1,J)*255. /FMAXU 
ELSE 
UCI, J)=C(UC1I,J)=FMINU)*255. /FMAXU 
ENDIF 
TOS CuCl) 
IF( IU. GT. 127) THEN 
IU=IU-256 
ENDIF 
BYTEU(J)-1U 
CONTINUE 
WRITE (3'I+1) (BYTEU(N) ,N=0,127) 
CONTINUE 
BOUSO0) 1—Cmelizg 
DO 850) J=Oniey 
IF(Y(1,J). LT. FMINY)THEN 
FMINY=Y(I,J) 
ENDIF 
IF(Y(1I,J). GT. FMAXY )THEN 
FMAXY=Y(I, I : 
ENDIF 
CONTINUE 
CONTINUE 
IF (FMINY. LT. 0) THEN 
FMAXY=FMAXY-FMINY 
ENDIF 
DO 900° IsOhie 
DO 950 J=0,127 
IF(FMINY. LT. 0) THEN 
Y(1,J)=Y(1,J)-FMINY 
YCI,J)=YC1,J)%255) / FMA 
ELSE 
Y(I,J}=C(Y(1,J)=FMINY)*255. /FMAXY 
ENDIF 
LY=NINT(YC1 
IF( IY. GT. 127) THEN 
IY=1Y-256 
ENDIF 
BYTEY(J)-I1Y 
CONTINUE 
WRITE (4'I+1) (BYTEY(N) ,N=0,127) 
CONTINUE 
CLOSE (UNIT=3) 
CLOSE (UNIT=4) 
STOP 
END 


SUBROUTINE GAUSS(IX,S,AM,V) 
AMP=0. 0 
DO 50 I=1,12 

RANDOM=RAN( IX) 


50 


50 


AMP=AMP+RANDOM 
CONTINUE 
V=(AMP-6. 0)*S+AM 
RETURN 
END 


SI 
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APPENDIX D. RLS SYSTEM IDENTIFICATION 


SYSTEM IDENTIFICATION OF AN ADAPTIVE RLS 2X2 FILTER VERSUS 
2X2 RLS FILTER VERSUS A A 2K2 FILTER WITH 
KNOWN COEFFICIENTS 


DIMENSION WC -1: 32,-1: 32) ,PC0: 3,0: 3) (AA 0: oe) 


KKKKEKVARTABLE DEFINITIONS **%¥%x%% 
A=NOISE AMPLITUDE 

S=NOISE VARIANCE 

AM=NOISE MEAN 

IX=SEED 

AA=FILTER COEFFICIENTS 

W=WHITE GAUSSIAN NOISE GENERATED BY SUBROUTINE 
D=DESIRED OUTPUT 

Y=FILTER OUTPUT 

E=ERROR(D - Y) 

MXM = ARRAY SIZE 

NXN = FILTER SIZE 

Q=WEIGHTING FACTOR 

P=INVERSE OF COVARIANCE MATRIX 


wededesereeee TNPUT INITIAL VALUES **%xtkvee eee 
M=32 
N=2 
A=1. 


Hou ual 
OOo © 


AM=0 
IX=65539 


Sessecceccee BUILD INITIAL ‘De MATRIX desevckdtcHeseadk 
DO 1 K=0,(N**2-1) 
DO 5 J=0,(N**2-1) 
IF(K. EQ. J) THEN 
P(K,J)=100. 


CONTINUE 
OPEN (UNIT=63,FILE='RLSID2D2 DATA' ,STATUS='OLD' ) 


weedeteeere CALCULATE INPUT MATRIX (WHITE GAUSSIAN NOISE) **%%#%* 
DO 10 K=-1,M-1 
DO 20 J=-1,M-1 
IF((K. LT. 0). OR. (J. LT. 0) ) THEN 
W(K,J)=0. 


By 


20 


moacaa 


50 
40 
15 
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20 
30 


40 


30 
60 


70 


80 


ELSE 
CALL GAUSS(IX,S,AM,V) 
WOK, J)=V*A 
ENDIF 
CONTINUE 
CONTINUE 


kik CALCULATE ERROR BETWEEN ADAPTIVE FILTER OUTPUT 

AND KNOWN FILTER OUTPUT. USE THIS ERROR TO UPDATE 

FILTER COEFFICIENTS AND THE 'P' MATRIX 9 kiiedededoederededeies 
DO 40 K=0,M-1 

DO 50 J=0,M-1 
L=K*M+J 
D=. 4*W(K,J)+. 6*W(K-1,J)-. 3*W(K, J-1)+. 2*W(K-1,J-1) 
Y=AA(0)*W(K,J)+AA(1)*W(K-1,J)+AA(2)*W(K,J-1)+ 
AA(3)*W(K-1,J-1) 

E=D-Y 

X(0)=W(K,J) 

X(1)=W(K-1,J) 

X(2)=W(K,J-1) 

X(3)=W(K-1,J-1) 

PAGER LS Cr, x,AA,N,E,Q) 
PRINT 15, L,E,AA(0),AA(1),AA(2) , AA(3) 
WRITE (63,15) L,E,AA(0) ,AA(1),AA(2) , AA(3) 

CONTINUE 
CONTINUE 
FORMAT (' ',2X,13,2X,F9.5,2X,F9.5,2X,F9.5,2X,F9.5,2X,F9. 5) 
STOP 
END 


Tekkewkeekeee = =TO COMPUTE THE GAIN MATRIX xxv 
SUBROUTINE RLS(CP,X,AA,N,E,Q) 
DIMENSION P(0: 3,0: 3),X(0: 3) ,AACO: 3), TEMY(0: 3,0: 3), 


X GAMA(0: 3) , TEMC(0: 3) 


DO 30 L=0,3 
TEM=0. 
DO 20 K=0,3 
TEM=TEM+X(K)*P(K,L) 
GAMA( L)=TEM 
GAM=0 
DO 40 K=0,3 
GAM=GAM+GAMA(K)*X(K) 
DO 60 L=0,3 
TEM=0. 
DO 50 K=0,3 
TEM=TEM+P(L,K)*X(K) 
TEMC( L)=TEM 
DO 70 L=0,3 
DO 70 K=0,3 
TEMY(L,K)=TEMC( L)*GAMA(K) 
DO 80 K=0,3 
DO 80 L=0,3 
P(K,L)=(P(K, L) -CTEMY(K,L) /GAM))/Q 
DO 100 K=0,3 
TEM=0, 


0 


20 
100 


110 


DOSI0) E=Or 2 
TEM=TEM+P(K,L)*X(L) 
TEMC(K)=TEM 
DO 110 K=0,3 
AA(K)=AA(CK)+TEMC(K)*E 
Q=. 99*Q+. O1 
RETURN 
END 
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APPENDIX E. RLS ALGORITHM IMPLEMENTED IN A NOISE 


CANCELER 


THIS IS A VAX/VMS FORTRAN PROGRAM THAT IMPLEMENTS A TWO 
DIMENSIONAL 2X2 ADAPTIVE RLS FILTER WITHIN A NOISE 
CANCELER 


INTEGER M,N,K,J 

BYTE B(0: 127) ,BYTEE(0: 127) ,BYTEU(0: 127) 

INTEGER*4 INTE(0: 127,0:127),INTU(O: 127,0: 127), IE, IU 
REAL*4 A,FMINE,FMAXE , FMINU, FMAXU 

REAL*4 AA(0: 3),E(0: 127,0: 127) ,U(0: 127,0: 127) 

REAL*4 P(0:3,0:3),X(0:3),IM(0: 127,0: 127) 

REAL*4 W(-1:127,-1:127),Y(0: 127,0: 127) 


HeRHKKVARTABLE DEFINITIONS *¥%r% 

A=NOISE AMPLITUDE 

=NOISE VARIANCE 

AM=NOISE MEAN 

AA=FILTER COEFFICIENTS 

IX=SEED 

IM=INPUT IMAGE 
W=WHITE GAUSSIAN NOISE GENERATED BY SUBROUTINE 
U=SIGNAL PLUS NOISE( IM + W) 

Y=FILTER OUTPUT 

E=ERROR(U - Y) 

=WEIGHTING FACTOR 

P=INVERSE OF COVARIANCE MATRIX 

FMINE , FMAXE , FMINU , FAMXU=PARAMETERS 

TO BE USED TO CONVERT DECIMAL DATA TO BYTE DATA 

MXM=ARRAY SIZE 

NXN=FILTER SIZE 


wee INITIAL VALUES xsevedaesedex 
M=128 

N=2 

Q=1. 

A=1. 

AAC 0)=0. 

AAC 1)=0. 

AAC 2)=0. 
AA(3)=0. 

S=40. 

AM=0. 
TX=65539 
FMINE=1. E+10 
FMAXE=-1. E+10 
FMINU=1. E+10 
FMAXU=-1.E+10 


aS 


202°) 


110 
100 


ONIN OrF WN 
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weieieievexeQPEN AN IMAGE FILE, CONVERT THE BYTE DATA INTO INTEGERS 
AND THEN PLACE THESE VALUES IN A MATRIX:iorddedex 
OPEN (UNIT=1, NAME ='HOUS1G.DAT', TYPE ='OLD', ACCESS = 
x 'DIRECT', RECORDSIZE=32, MAXREC=128) 
DO 100 K=0,127 
READ(1'K+1) (B(L),L=0,127) 
DO 110 J=0,127 
IF(B( J). LT. 0) THEN 
IM(K,J)=B(J)+256 
ELSE 
IM(K,J)=B(J) 
ENDIF 
CONTINUE 
CONTINUE 
CLOSE (UNIT=1) 


teokieieke BUILD INITIAL 'P' MATRIX xvvedercdedex 
DO 1 K=0,(N**2-1) 
DO 5 J=0,(N**2-1) 
IF CKeEOs J ) THEN 
P(K,J)=100. 
ELSE 
P(K,J)=0. 0 
ENDIF 
CONTINUE 
CONTINUE 


‘eke ADD WHITE GAUSSIAN NOISE TO THE IMAGE AND SET THE 
VALUES OUTSIDE THE IMAGE TO ZERO*ocddex 
DO 10 K=-1,M-1 
DO 20 J=-1,M-1 
IF ((K. LT. 0). OR. (J. LT. 0))THEN 
W(K,J)=0. 
IM(K,J)=0. 
ELSE 
CALL GAUSS(IX,S,AM,V) 
W(K, J)=V*A 
ENDIF 
UC(K, J)=IM(K, J)4+W(K,J) 
CONTINUE 
CONTINUE 


week KE CALCULATE ERROR BETWEEN ADAPTIVE FILTER OUTPUT 
AND KNOWN FILTER OUTPUT. USE THIS ERROR TO UPDATE 
FILTER COEFFICIENTS AND THE 'P' MATRIX 9X tevedevedevededededeiede 
DO 40 K=0,M-1 
DO 50 J=0,M-1 
L=K*M+J 
Y=AA(0)*W(K,J)+AA(1)*W(K-1,J)+AA(2)*W(K,J-1)+ 
* AA(3)*W(K-1,J-1) 
E(K,J)=U(K,J)-Y(K,J) 
EE=E(K,J) 
X(0)=W(K,J) 
X(1)=W(K-1,J) 
X(2)=W(K,J-1) 
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C 
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>> 0 
500 


650 


600 


750 
700 


MCS )=WCK-1,J-1) 
CALL RLS(P,X,AA,N,EE,Q) 
CONTINUE 
CONTINUE 


kere CHANGE SIGNAL PLUS NOISE AND ERROR OUTPUT 
INTO BYTE DATA AND WRITE TO A FILE***%*% 
OPEN (UNIT=2 ,NAME='ERROR. DAT’ , TYPE='NEW’ , ACCESS= 
‘DIRECT’ , RECORDSIZE=32 , MAXREC=128) 
OPEN (UNIT=3 , NAME='SIGNOISE. DAT' , TYPE='NEW' , ACCESS= 
"DIRECT ,RECORDSIZE=32 , MAXREC=128) 
DO 500 I=0,127 
DO 550 J=0,127 
IF(E(1,J). LT. FMINE) THEN 
FMINE=E(1I,J) 
ENDIF 
IF(E(1I,J). GT. FMAXE)THEN 
FMAXE=E(I,J) 
ENDIF 
CONTINUE 
CONTINUE 
IF (FMINE. LT. 0) THEN 
FMAXE=FMAXE -FMINE 
ENDIF 
DO 600 I=0,127 
Bowos0 J=0,127 
IF(FMINE. LT. 0)THEN 
E(I,J)=E(1,J)-FMINE 
E(1,J)=E(1,J)*255. /FMAXE 
ELSE 
E(I,J)=(E(1,J)=FMINE)*255. /FMAXE 
ENDIF 
IE=NINT(E(1I,J)) 
IFC IE. GT. 127) THEN 
IE=IE-256 
ENDIF 
BYTEE(J)-IE 
CONTINUE 
WRITE (2'I+1) (BYTEE(N) ,N=0,127) 
CONTINUE 
Dor 00 I=0,127 
BO 750 J-0, 1127 
IF(U(I,J). LT. FMAINU)THEN 
FMINU=U(I,J) 
ENDIF 
IF(UCI,J). GT. FMAXU)THEN 
FMAXU=U(I1,J) 
ENDIF 
CONTINUE 
CONTINUE 
IF (FMINU. LT. 0) THEN 
FMAXU=FMAXU-FMINU 
ENDIF 
DOe400 1=0,127 
Bo 450 J=0,127 
IF(FMINU. LT. 0)THEN 
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450 


400 


20 
30 


40 


510, 
60 


70 


80 


90 
100 


110 


UCI, J)=UC I ,J)-FMINU 
UCI, J)=UC1I,J)*255. /FMAXU 
ELSE 
UC I,J)=CUCI,J)=FMINU)*255. /FMAXU 
ENDIF 
IU=NINT(U(CI,J)) 
IF(IU.GT. 127) THEN 
IU=IU-256 
ENDIF 
BYTEU( J) -1U 
CONTINUE 
WRITE (3'I+1) (BYTEUCN) ,N=0,127) 
CONTINUE 
CLOSE (UNIT=2) 
CLOSE (UNIT=3) 
STOP 
END 


kite = =TO COMPUTE THE GAIN MATRIX *xX%XKK KKK 
SUBROUTINE RLS(P,X,AA,N,E,Q) 
DIMENSION P(0: 3,0: 3),X(0:3),AAC0: 3), TEMY(0: 3,0: 3), 
GAMA(C0: 3) , TEMCC0: 3) 
DO 30 L=0,3 
TEM=0. 
DO enk=o> 
TEM=TEM+X(K)*P(K, L) 
GAMA( L)=TEM 
GAM=Q 
DO 40 K=0,3 
GAM=GAM+GAMA(C K)*X(K) 
DO 60 L=0,3 
TEM=0. 
DO 50 K=0,3 
TEM=TEM+P(L,K)*X(K) 
TEMC( L)=TEM 
DO 70er—ons 
DO 70 K=0,3 
TEMY(L,K)=TEMC( L)*GAMA(K) 
DO 80 K=0,3 
DO 80el—oes 
PCR SUE CK rhe Gan) 1/ © 
DO 100 K=0,3 
TEM=0. 
DO 90 L=0,3 
TEM=TEM+P(K,L)*X(CL) 
TEMC(K)=TEM 
BO Women 
AACK)=AACK)+TEMCCK)*E 
Q=. 99*Q+. 01 
RETURN 
END 
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APPENDIX F. RLS ALGORITHM IMPLEMENTED IN AN ADAPTIVE 
LINE ENHANCER 


Gla Ca C2 


©y C2°C) C2). C2 CD Clearances C2)Ca C2 CIC IC) Cac) C2 CD 


THIS IS A VAX/VMS FORTRAN PROGRAM THAT IMPLEMENTS A TWO APPO0O0010 
DIMENSIONAL 2X2 ADAPTIVE RLS FILTER WITHIN AN ADAPTIVE APP00020 
LINE ENHANCER APP00030 
APP00040 

INTEGER M,N,K,J APPOO050 
Bynes 2/ ),BYTEUCO: 127) ,BYTEYCO: 127) APP00060 
Piece RS4 INTY(0: 127,0: 127), INTUCO: 127,0: 127) ,IU,IY APP00070 
REAL*4 A,FMINU, FMAXU, FMINY , FMAXY APP00080 
REAL*4 AA(0:3),E(0: 127,0: 127) ,U(0: 127,0: 127) APP00090 
Reese PCO: 3,0:3),X(0:3),1M(C0: 127,0: 127) APP00100 
REAL*4 W(-1:127,-1:127),Y(0: 127,0: 127) APP00110 
APEUOIZY) 

HIKRKAKVARTABLE DEFINITIONS ¥ Xt APP00130 
A=NOISE AMPLITUDE APP00140 
S=NOISE VARIANCE APE OOS 0 
AM=NOISE MEAN APP00160 
AA=FILTER COEFFICIENTS APP00170 
IX=SEED APP00180 
IM=INPUT IMAGE APP00190 
=WHITE GAUSSIAN NOISE GENERATED BY SUBROUTINE APP00200 
U=SIGNAL PLUS NOISECIM + W) APP00210 
Y=FILTER OUTPUT APP00220 
E=ERROR(U - Y) APP00230 
WW=DELAYED VERSION OF SIGNAL PLUS NOISE(U) APP00240 
=WEIGHTING FACTOR APP00250 
P=INVARSE OF COVARIANCE MATRIX APP00260 
FMINU , FAMXU,FMINY , FMAXY=PARAMETERS APP00270 
TO BE USED TO CONVERT DECIMAL DATA TO BYTE DATA APP00280 
MXM=ARRAY SIZE APP00290 
NXN=FILTER SIZE APP00300 
APP00310 

eee TNITILAL VALUE Stvccedecededed APP00320 
M=128 APP00330 
N=2 APP00340 
Q=1. APP00350 
A=1. APP00360 
AA(0)=0. APP00370 
AA(1)=0. APP00380 
AA(2)=0. APP00390 
AA(3)=0. APP00400 
S=40. APP00410 
AM=0O. APP00420 
IX=65539 APP00430 
FMINU=1. E+10 APP00440 
FMAXU=-1. E+10 APP00450 
FMINY=1. E+10 APP00460 
FMAXY=-1. E+10 APP00470 
APP00480 


ay 


110 
100 


Qa ae MN 


46 
45 


eek OPEN AN IMAGE FILE, CONVERT THE BYTE DATA INTO INTEGERS 
AND THEN PLACE THESE VALUES IN A MATRIX*tetededevex 
OPEN (UNIT=1, NAME ='HOUS1G. DAT', TYPE ='OLD', ACCESS = 


* "DIRECT', RECORDSIZE=32, MAXREC=128) 


DO 100 K=0,127 
READ(1'K+1) (BCL) ,L=0,127) 
DO 110 J=0,127 
IF(B( J). LT. 0) THEN 
IM(K, J)=B(J)+256 
ELSE 
IM(K,J)=B(J) 
ENDIF 
CONTINUE 
CONTINUE 
CLOSE (UNIT=1) 


tomer BULLD INITIAL 'P' MATRIX xk 
DO 1 K=0,(N**2-1) 
DO 5 J=0,(N**2-1) 
IF(K. EQ. J) THEN 
P(K,J)=100. 
ELSE 
P(K,J)=0. 0 


CONTINUE 
CONTINUE 


dexvevexeteeADD WHITE GAUSSIAN NOISE TO THE IMAGE AND SET THE 
VALUES OUTSIDE THE IMAGE TO ZERO*iienke 
DO 10 K=-2,M-1 
DO 20 J=-2,M-1 
IF ((K. LT. 0). OR. (J. LT. 0))THEN 
W(K,J)=0. 
IM(K,J)=0. 
ELSE 
CALL GAUSS(IX,S,AM,V) 
W(K, J)=V*A 
ENDIF 
U(K, J)=IM(K, J)+W(K, J) 
CONTINUE 
CONTINUE 


week COMPUTE THE DELAYED VERSION OF THE SIGNAL 
PLUS NOISE MATRIX(U) 7% 
DO 45 K=-1,127 
DO 46 J=-1,127 
IF((M. LT. -1). OR. (N. LT. -1) THEN 
WW(K,J)=0. 
ELSE 
WW(K,J)=U(K-1,J) 
ENDIF 
CONTINUE 
CONTINUE 


veered CALCULATE ERROR BETWEEN ADAPTIVE FILTER OUTPUT 
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APP00490 
APP00500 
APP00510 
APP00520 
APP00530 
APP00540 
APP00550 
APP00560 
APP00570 
APP00580 
APP00590 
APP00600 
APP00610 
APP00620 
APP00630 
APP00640 
APP00650 
APP00660 
APP00670 
APP00680 
APP00690 
APP00700 
APP00710 
APPO0720 
APP00730 
APP00740 
APP00750 
APP00760 
APP00770 
APP00780 
APP00790 
APP00800 
APP00810 
APP00820 
APP00830 
APP00840 
APP00850 
APP00860 
APP00870 
APP00880 
APP00890 
APP00900 
APP00910 
APP00920 
APP00930 
APP00940 
APP00950 
APP00960 
APP00970 
APP00980 
APP00990 
APPO1000 
APP01010 
APP01020 
APPO01030 


50 
40 


750 
700 


450 


400 


AND KNOWN FILTER OUTPUT. USE THIS ERROR TO UPDATE 
FILTER COEFFICIENTS AND THE 'P' MATRIX sievedededrieedeiriesex 
DO 40 K=0,M-1 
DO 50 J=0,M-1 
Y=AA(0)*WW(K,J)+AA(1)*WW(K-1,J)+AA(2)*WW(K,J-1)+ 


% AA(3)*WW(K-1,J-1) 
E(K,J)=U(K,J)-Y(K,J) 
EE=E(K,J) 


X(0)=WW(K, J) 
X(1)=WW(K-1,J) 
X(2)=WW(K,J-1) 
X(3)=WW(K-1,J-1) 
GAME RLS(P,X,AA,N,EE,Q) 
CONTINUE 
CONTINUE 


HIKARI CHANGE SIGNAL PLUS NOISE AND FILTER OUTPUT 
INTO BYTE DATA AND WRITE TO A FILE*#%%% 
OPEN (UNIT=3 , NAME=' SIGNOISE. DAT’ , TYPE='NEW' , ACCESS= 
* "DIRECT , RECORDSIZE=32 , MAXREC=128) 
OPEN (UNIT=4 , NAME= "FILTERED. DAT’ , TYPE='NEW' , ACCESS= 
* "DIRECT , RECORDSIZE=32 , MAXREC=128) 
DO 700 I=0,127 
POs750 J=0, 127 
IF(UC I,J). LT. FMINU) THEN 
FMINU=U(I, J) 
ENDIF 
IF(U( I,J). GT. FMAXU) THEN 
FMAXU=U(I,J) 
ENDIF 
CONTINUE 
CONTINUE 
IF (FMINU. LT. 0) THEN 
FMAXU=FMAXU-FMINU 
ENDIF 
moecoO I=0,127 
BONS O10), 127 
IF(FMINU. LT. 0)THEN 
U(1,J)=U( I,J) -FMINU 
UCI, J)=UC1I ,J)*255. /FMAXU 
ELSE 
UCI, J)=(UC1I,J)=FMINU)*255. /FMAXU 
ENDIF 
IU=NINT(U(I,J)) 
IF( IU. GT. 127) THEN 
IU=IU-256 
ENDIF 
BYTEU( J)-IU 
CONTINUE 
WRITE (3'I+1) (BYTEU(N) ,N=0,127) 
CONTINUE 
DO 800 I=0,127 
BOwcs 0) 3=0),127 
IF(Y( I,J). LT. FMINY) THEN 
FMINY=Y(I,J) 
ENDIF 
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APP01040 
APP0O1050 
APP01060 
APP01070 
APP01080 
APP01090 
APP01100 
APP01110 
APP01120 
APP01130 
APP01140 
APP0O1150 
APP0O1160 
APP01170 
APP01180 
APP0O1190 
APP01200 
APP01210 
APP01220 
APP01230 
APP01240 
APP01250 
APP01260 
APP01270 
APP01280 
APP01290 
APP01300 
APP01310 
APP01320 
APP01330 
APP01340 
APP01350 
APP01360 
APP01370 
APP01380 
APP01390 
APP01400 
APP01410 
APP01420 
APP01430 
APP01440 
APP01450 
APP01460 
APP01470 
APP01480 
APP01490 
APP01500 
APPO1510 
REPO 15> Z0 
AEEOIS 50 
APP01540 
APP01550 
APP01560 
APP01570 
APP0O1580 
APP01590 


850 
800 


20 


900 


20 
30 


40 


50 
60 


70 


80 


IF(Y(1,J). GT. FMAXY)THEN 
FMAXY=Y(I,J) 
ENDIF 
CONTINUE 
CONTINUE 
IF (FMINY. LT. 0) THEN 
FMAXY=FMAXY -FMINY 
ENDIF 
DO 900 I=0,127 
DO 950 J=0,127 
IF(FMINY. LT. 0) THEN 
Y(1I,J)=Y(1,J)-FMINY 
Y(I,J)=Y(1,J)*255. /FMAXY 
ELSE 
Y(I,J)=(Y(1,J)=FMINY)*255. /FMAXY 
ENDIF 
IY=NINT(Y(I,J)) 
IF(IY.GT.127) THEN 
IY=IY-256 
ENDIF 
BYTEY(J)-IY 
CONTINUE 
WRITE (4'I+1) (BYTEY(N) ,N=0,127) 
CONTINUE | 
CLOSE (UNIT=3) 
CLOSE (UNIT=4) 
STOP 
END 


tiskkieeke TO COMPUTE THE GAIN MATRIX *vevededededededd 
SUBROUTINE RLS(P,X,AA,N,E,Q) 
DIMENSION P(0: 3,0: 3),X(C0: 3) ,AAC O03), TEMY(0: 30 
GAMA(0: 3) , TEMC(0: 3) 
DO 30 L=0,3 
TEM=0. 
DO 20 K=0,3 
TEM=TEM+X(K)*PCK,L) 
GAMA( L)=TEM 
GAM=Q 
DO 40 K=0,3 
GAM=GAM+GAMA( K)*X(K) 
DO 60 L=0,3 
TEM=0. 
DO 50 K=0,3 
TEM=TEM+P(L, K)*X(K) 
TEMC( L)=TEM 
DO 70 L=0,3 
DO 70 K=0,3 
TEMY(L, K)=TEMC( L)*GAMA(K) 
DO 80 K=0,3 
DO 80 L=0,3 
P(K,L)=(P(K, L)-(TEMY(K,L)/GAM) )/Q 
DO 100 K=0,3 
TEM=0. 
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APP01600 
APP0O1610 
APP01620 
APP01630 
APP01640 
APP01650 
APP01660 
APP01670 
APP01680 
APP01690 
APP01700 
APP0O1710 
APP01720 
APP01730 
APP01740 
APP01750 
APP01760 
APP01770 
APP01780 
APP01790 
APP01800 
APP01810 
APP01820 
APP01830 
APP01840 
APP01850 
APP01860 
APP01870 
APP01880 
APP01890 
APP01900 
APP01910 
APP01920 
APP01930 
APP01940 
APP01950 
APP01960 
APP01970 
APP01980 
APP01990 
APP02000 
APP02010 
APP02020 
APP02030 
APP02040 
APP02050 
APP02060 
APP02070 
APP02080 
APP02090 
APP02100 
APP02110 
APP02120 
APP02130 


90 
100 


110 


DO 90 L=0,3 
TEM=TEM+P(K,L)*X(L) 
TEMC( K)=TEM 
DO 110 K=0,3 
AA(K)=AA(K)+TEMC(K)*E 
Q=. 99*Q+. 01 
RETURN 
END 
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APP02140 
APP02150 
APP02160 
italy 
APP02180 
AREOZ 190 
APP02200 
Ar 2210 
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